1 Preparations

Load libraries and general settings:

## load libraries
library(tidyverse)
library(lme4)
library(lmerTest)
library(optimx)
library(emmeans)
library(psych)
library(ggpubr)
library(here)
emm_options(lmer.df = "asymptotic") # as pre-registered
             
two_colors <- c("#D55E00", "#56B4E9") # facilitation and interference
dark_colors <- c("#800D00", "#165F83", "gray45") # darker (facilitation, interference, null)
subt_shape <- 20 # diamond
regr_shape <- 18 # round
equi_boundary <- .15

# APA theme for figures
theme_set(papaja::theme_apa(base_size = 12, base_family = "Helvetica", box = FALSE))
theme_update(strip.placement = "outside")

1.1 Load data

# To run the code in this chunk, you need to put all raw data files into data/.

# # read all raw data, remove subjid and re-save as a new file
# df_raw_tmp <- list.files("data", pattern = "*.csv", full.names = TRUE) %>%
#   map_dfr(read_csv, show_col_types = FALSE, .id="Subject") %>% # overwrite Subject ID from Prolific
#   mutate(Subject = sprintf("subj%03d", as.integer(Subject))) # add new Subject
#
# write_rds(df_raw_tmp, here("data", "2faces_pro_raw480.rds"))
# write_csv(df_raw_tmp, here("data", "2faces_pro_raw480.csv"))

# connections between prolific subject ID and ID in this study
# df_subjects <- list.files("data", pattern = "*.csv", full.names = TRUE) %>%
#   map_dfr(read_csv, show_col_types = FALSE, .id="SubjectThis") %>% # overwrite Subject ID from Prolific
#   mutate(SubjectThis = sprintf("subj%03d", as.integer(SubjectThis))) %>%  # add new Subject
#   select(ProlificID = Subject, Subject = SubjectThis) %>% 
#   distinct()
# 
# write_rds(df_subjects, here("data", "2faces_pro_subjectID.rds"))
# write_csv(df_subjects, here("data", "2faces_pro_subjectID.csv"))

The structure of the raw data df_raw:

# load data
df_raw <- readRDS(here("data", "2faces_pro_raw480.rds")) %>% 
  filter(trial_frame=="test_face",
         Section=="main") %>% # exclude practice trials
  select(Subject, CBcode, Task_name:SameDifferent, 
         Correct_ans_posi, Correct_response, response, 
         Correct, RT, StimGroup, StudyFace, TestFace) %>% 
  transmute(Subject = as_factor(Subject),
            Task = as_factor(Task_name), # Task_order
            Trial_num,
            PW = as_factor(PW), 
            Feature = as_factor(Feature),
            Cue = as_factor(Cue),
            Congruency = as_factor(str_sub(Congruency, 1, 3)),
            Alignment = as_factor(str_sub(Alignment, 1, 3)),
            SD = as_factor(str_sub(SameDifferent, 1, 4)),
            SD = fct_relevel(SD, "same", "diff"),
            CorrectAnswer = Correct_response, 
            Response = response,
            Correct,
            RT,
            StimGroup = if_else(Task_name=="PW", 
                                str_remove(basename(StudyFace),".png"), 
                                StimGroup),
            StudyFace, 
            TestFace) 

str(df_raw)
## tibble [337,920 × 16] (S3: tbl_df/tbl/data.frame)
##  $ Subject      : Factor w/ 480 levels "subj001","subj002",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ Task         : Factor w/ 3 levels "SCF","PW","CCF": 1 1 1 1 1 1 1 1 1 1 ...
##  $ Trial_num    : num [1:337920] 1 2 3 4 5 6 7 8 9 10 ...
##  $ PW           : Factor w/ 2 levels "whole","part": NA NA NA NA NA NA NA NA NA NA ...
##  $ Feature      : Factor w/ 3 levels "mouth","nose",..: NA NA NA NA NA NA NA NA NA NA ...
##  $ Cue          : Factor w/ 2 levels "top","bot": 1 1 1 1 1 1 1 1 1 1 ...
##  $ Congruency   : Factor w/ 3 levels "con","inc","iso": 1 1 1 2 2 1 2 1 2 2 ...
##  $ Alignment    : Factor w/ 3 levels "ali","mis","iso": 1 1 2 2 2 1 2 2 2 1 ...
##  $ SD           : Factor w/ 2 levels "same","diff": 2 2 2 1 1 2 1 2 1 1 ...
##  $ CorrectAnswer: chr [1:337920] "1" "1" "1" "2" ...
##  $ Response     : chr [1:337920] "2" "2" "2" "1" ...
##  $ Correct      : logi [1:337920] FALSE FALSE FALSE FALSE FALSE FALSE ...
##  $ RT           : num [1:337920] 1117 810 1793 719 721 ...
##  $ StimGroup    : chr [1:337920] "M4" "F5" "F3" "F4" ...
##  $ StudyFace    : chr [1:337920] "M4197_M4136_ali_top.png" "F5068_F5078_ali_top.png" "F3047_F3079_mis_top.png" "F4064_F4018_mis_top.png" ...
##  $ TestFace     : chr [1:337920] "M4134_M4188_ali_top.png" "F5076_F5033_ali_top.png" "F3112_F3062_mis_top.png" "F4064_F4040_mis_top.png" ...

1.2 Tidy data

The information on inclusion and exclusion criteria can be found in the online pre-registration.

1.2.1 Exclude participants

We will first exclude participants who had 10% or more trials in any of the three tasks with either very fast or very slow responses (<200ms or >5000ms).

df_ex_subj <- df_raw %>% 
  mutate(isout = RT < 200 | RT > 5000) %>%  # identify outlier trials
  group_by(Subject, Task) %>% 
  summarize(count_all = n(),
            count_out = sum(isout),
            .groups = "drop") %>% 
  mutate(ratio = count_out/count_all) %>% 
  filter(ratio >= .1) # who had 10% or more trials 

(subjlist_ex <- unique(df_ex_subj$Subject))
##  [1] subj034 subj045 subj074 subj094 subj123 subj143 subj150 subj154 subj171 subj188 subj198 subj204 subj241 subj244 subj255 subj262 subj270 subj288 subj289 subj314 subj335 subj383 subj440 subj459 subj461
## 480 Levels: subj001 subj002 subj003 subj004 subj005 subj006 subj007 subj008 subj009 subj010 subj011 subj012 subj013 subj014 subj015 subj016 subj017 subj018 subj019 subj020 subj021 subj022 subj023 subj024 subj025 subj026 subj027 subj028 subj029 subj030 subj031 subj032 subj033 subj034 subj035 subj036 subj037 subj038 subj039 subj040 subj041 subj042 subj043 subj044 subj045 subj046 subj047 subj048 subj049 subj050 subj051 subj052 subj053 subj054 subj055 subj056 subj057 subj058 subj059 subj060 subj061 subj062 subj063 subj064 subj065 subj066 subj067 subj068 subj069 subj070 subj071 subj072 subj073 subj074 subj075 subj076 subj077 subj078 subj079 subj080 subj081 subj082 subj083 subj084 subj085 subj086 subj087 subj088 subj089 subj090 subj091 subj092 subj093 subj094 subj095 subj096 subj097 subj098 subj099 subj100 subj101 subj102 subj103 subj104 subj105 subj106 subj107 subj108 subj109 subj110 subj111 subj112 subj113 subj114 subj115 subj116 subj117 subj118 subj119 subj120 subj121 subj122 subj123 subj124 subj125 subj126 subj127 subj128 subj129 subj130 subj131 subj132 subj133 subj134 subj135 subj136 subj137 subj138 subj139 subj140 subj141 subj142 subj143 subj144 subj145 subj146 subj147 subj148 subj149 subj150 subj151 subj152 subj153 subj154 subj155 subj156 subj157 subj158 subj159 subj160 subj161 subj162 subj163 subj164 subj165 subj166 subj167 subj168 subj169 subj170 subj171 subj172 subj173 subj174 subj175 subj176 subj177 subj178 subj179 subj180 subj181 subj182 subj183 subj184 ... subj480

In summary, 25 participants were excluded, and therefore, there are 455 remaining (valid) participants.

1.2.2 Exclude trials

Second, we exclude outlier trials in the remaining participants. We calculate the standard deviations (SDs) of response times (RTs) for each condition and each participant in each task separately (e.g., for part-whole task, such a condition would be a whole-face condition with a specific target feature such as the eyes. For the standard composite task, such a condition would be same-aligned trials). Trials with responses made <200ms or > 3SDs of the condition will be excluded.

On average, the percentages of trials removed for each task:

df_raw <- df_raw %>% 
  filter(!Subject %in% subjlist_ex) %>% # exclude outlier subjects
  group_by(Subject, Task, PW, Feature, Cue, Congruency, Alignment, SD) %>% 
  summarize(Trial_num, CorrectAnswer, Response, Correct, RT, 
            StimGroup, StudyFace, TestFace,
            Z_rt = scale(RT), # calculate Z value
            .groups = "drop") %>% 
  mutate(isout = Z_rt < -3 | Z_rt > 3 | RT < 200) 

df_raw %>% 
  group_by(Subject, Task) %>% 
  summarize(count = n(),
            count_out = sum(isout),
            .groups = "drop") %>% 
  mutate(ratio = count_out/count) %>% 
  group_by(Task) %>% 
  summarize(count = mean(count),
            avg_outlier_N = mean(count_out),
            `avg_outlier_ratio (%)` = mean(ratio)*100) 

1.3 Save data for each task separately

The structure of the clean data (df_tidy):

df_tidy <- df_raw %>% 
  filter(!isout) %>% 
  select(-c(Z_rt, isout))

str(df_tidy)
## tibble [314,222 × 16] (S3: tbl_df/tbl/data.frame)
##  $ Subject      : Factor w/ 480 levels "subj001","subj002",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ Task         : Factor w/ 3 levels "SCF","PW","CCF": 1 1 1 1 1 1 1 1 1 1 ...
##  $ PW           : Factor w/ 2 levels "whole","part": NA NA NA NA NA NA NA NA NA NA ...
##  $ Feature      : Factor w/ 3 levels "mouth","nose",..: NA NA NA NA NA NA NA NA NA NA ...
##  $ Cue          : Factor w/ 2 levels "top","bot": 1 1 1 1 1 1 1 1 1 1 ...
##  $ Congruency   : Factor w/ 3 levels "con","inc","iso": 1 1 1 1 1 1 1 1 1 1 ...
##  $ Alignment    : Factor w/ 3 levels "ali","mis","iso": 1 1 1 1 1 1 1 1 1 1 ...
##  $ SD           : Factor w/ 2 levels "same","diff": 2 2 2 2 2 2 2 2 2 2 ...
##  $ Trial_num    : num [1:314222] 1 2 6 13 23 32 35 36 38 43 ...
##  $ CorrectAnswer: chr [1:314222] "1" "1" "1" "1" ...
##  $ Response     : chr [1:314222] "2" "2" "2" "2" ...
##  $ Correct      : logi [1:314222] FALSE FALSE FALSE FALSE TRUE TRUE ...
##  $ RT           : num [1:314222] 1117 810 786 718 593 ...
##  $ StimGroup    : chr [1:314222] "M4" "F5" "M5" "M4" ...
##  $ StudyFace    : chr [1:314222] "M4197_M4136_ali_top.png" "F5068_F5078_ali_top.png" "M5178_M5146_ali_top.png" "M4136_M4134_ali_top.png" ...
##  $ TestFace     : chr [1:314222] "M4134_M4188_ali_top.png" "F5076_F5033_ali_top.png" "M5152_M5154_ali_top.png" "M4188_M4197_ali_top.png" ...
# custom function to set contrasts for factors
set_sdif <- function(.data, colstr, n=2){
  # .data dataframe
  # colstr the string of the column name
  
  contrasts(.data[[colstr]]) <- MASS::contr.sdif(n)
  
  return(.data)
}

1.3.1 Part-whole task

The structure of the clean data for part-whole task (df_pw):

df_pw <- df_tidy %>% 
  filter(Task == "PW",
         Feature != "nose") %>% # do not include Nose trials in GLMM
  select(Subject, PW, Feature, Correct, RT, StimGroup, StudyFace, TestFace) %>% 
  mutate(PW = fct_relevel(PW, "whole", "part"),
         Feature = fct_drop(Feature),
         PW_C = if_else(PW=="whole", -.5, if_else(PW=="part", .5, NaN)),
         Feature_C = if_else(Feature=="mouth", -.5, 
                             if_else(Feature=="eyes", .5, NaN)),
         PW_Feature = PW_C * Feature_C) %>% 
  set_sdif("PW") %>% 
  set_sdif("Feature")

# saveRDS(df_pw, file=here("jubail", "df_pw.rds"))

str(df_pw)
## tibble [42,880 × 11] (S3: tbl_df/tbl/data.frame)
##  $ Subject   : Factor w/ 480 levels "subj001","subj002",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ PW        : Factor w/ 2 levels "whole","part": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "whole" "part"
##   .. .. ..$ : chr "2-1"
##  $ Feature   : Factor w/ 2 levels "mouth","eyes": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "mouth" "eyes"
##   .. .. ..$ : chr "2-1"
##  $ Correct   : logi [1:42880] TRUE FALSE FALSE TRUE TRUE TRUE ...
##  $ RT        : num [1:42880] 2876 3001 1249 1439 1696 ...
##  $ StimGroup : chr [1:42880] "cau_m3" "cau_m5" "cau_m4" "cau_f5" ...
##  $ StudyFace : chr [1:42880] "stim/pw/main_stim/cau_m3.png" "stim/pw/main_stim/cau_m5.png" "stim/pw/main_stim/cau_m4.png" "stim/pw/main_stim/cau_f5.png" ...
##  $ TestFace  : chr [1:42880] "stim/pw/main_stim/cau_m3.mouth.whole.left.png" "stim/pw/main_stim/cau_m5.mouth.whole.left.png" "stim/pw/main_stim/cau_m4.mouth.whole.left.png" "stim/pw/main_stim/cau_f5.mouth.whole.right.png" ...
##  $ PW_C      : num [1:42880] -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 ...
##  $ Feature_C : num [1:42880] -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 ...
##  $ PW_Feature: num [1:42880] 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 ...

1.3.2 Standard composite face task

The structure of the clean data for standard composite face task (df_scf):

df_scf <- df_tidy %>% 
  filter(Task == "SCF",
         SD == "same") %>% # only use trials when the top is the same
  select(Subject, Alignment, Correct, RT, Trial_num, StimGroup) %>% 
  mutate(Alignment = fct_drop(Alignment),
         Ali_C = if_else(Alignment=="ali", -.5, 
                         if_else(Alignment=="mis", .5, NaN))) %>% 
  set_sdif("Alignment")

# saveRDS(df_scf, file=here("jubail", "df_scf.rds"))

str(df_scf)
## tibble [35,648 × 7] (S3: tbl_df/tbl/data.frame)
##  $ Subject  : Factor w/ 480 levels "subj001","subj002",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ Alignment: Factor w/ 2 levels "ali","mis": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "ali" "mis"
##   .. .. ..$ : chr "2-1"
##  $ Correct  : logi [1:35648] FALSE TRUE TRUE TRUE TRUE TRUE ...
##  $ RT       : num [1:35648] 749 528 817 798 1729 ...
##  $ Trial_num: num [1:35648] 10 20 22 25 27 28 31 39 41 42 ...
##  $ StimGroup: chr [1:35648] "F5" "F4" "F2" "F1" ...
##  $ Ali_C    : num [1:35648] -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 ...

1.3.3 Complete composite face task

The structure of the clean data for complete composite face task (df_ccf):

# including isolated condition
df_ccf_iso <- df_tidy %>% 
  filter(Task == "CCF") %>% # only use trials when the top is the same
  mutate(isSame = if_else( # whether participants responded "same" (signal)
    ((SD=="same")&Correct) | ((SD=="diff")&!Correct), 1,
    if_else(((SD=="same")&!Correct) | ((SD=="diff")&Correct), 0, NaN))) %>% 
  select(Subject, Cue, Congruency, Alignment, SD, Correct, 
         isSame, RT, Trial_num, StimGroup) 

# df_ccf for composite effect
df_ccf <- df_ccf_iso %>% 
  filter(Alignment != "iso") %>% 
  mutate(Alignment = fct_drop(Alignment),
         Congruency = fct_drop(Congruency)) %>% 
  set_sdif("Cue") %>% 
  set_sdif("Congruency") %>% 
  set_sdif("Alignment") %>% 
  set_sdif("SD") %>% 
  mutate(Con_C = if_else(Congruency=="con", -.5,
                         if_else(Congruency=="inc", .5, NaN)),
         Ali_C = if_else(Alignment=="ali", -.5,
                         if_else(Alignment=="mis", .5, NaN)),
         SD_C = if_else(SD=="same", -.5, if_else(SD=="diff", .5, NaN)),
         Con_Ali = Con_C * Ali_C,
         Con_SD = Con_C * SD_C,
         Ali_SD = Ali_C * SD_C,
         Con_Ali_SD = Con_Ali * SD_C)

# saveRDS(df_ccf, file=here("jubail", "df_ccf.rds"))

str(df_ccf)
## tibble [143,002 × 17] (S3: tbl_df/tbl/data.frame)
##  $ Subject   : Factor w/ 480 levels "subj001","subj002",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ Cue       : Factor w/ 2 levels "top","bot": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "top" "bot"
##   .. .. ..$ : chr "2-1"
##  $ Congruency: Factor w/ 2 levels "con","inc": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "con" "inc"
##   .. .. ..$ : chr "2-1"
##  $ Alignment : Factor w/ 2 levels "ali","mis": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "ali" "mis"
##   .. .. ..$ : chr "2-1"
##  $ SD        : Factor w/ 2 levels "same","diff": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "same" "diff"
##   .. .. ..$ : chr "2-1"
##  $ Correct   : logi [1:143002] TRUE TRUE TRUE TRUE TRUE TRUE ...
##  $ isSame    : num [1:143002] 1 1 1 1 1 1 1 1 1 0 ...
##  $ RT        : num [1:143002] 651 1399 986 819 1362 ...
##  $ Trial_num : num [1:143002] 56 70 76 90 113 114 122 146 147 151 ...
##  $ StimGroup : chr [1:143002] "F1" "M5" "F2" "M4" ...
##  $ Con_C     : num [1:143002] -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 ...
##  $ Ali_C     : num [1:143002] -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 ...
##  $ SD_C      : num [1:143002] -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 ...
##  $ Con_Ali   : num [1:143002] 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 ...
##  $ Con_SD    : num [1:143002] 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 ...
##  $ Ali_SD    : num [1:143002] 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 ...
##  $ Con_Ali_SD: num [1:143002] -0.125 -0.125 -0.125 -0.125 -0.125 -0.125 -0.125 -0.125 -0.125 -0.125 ...

The structure of the clean data for analysis involving isolated conditions for calculating facilitation and interference (df_fi):

# facilitation and interference

# df for facilitation and interference
df_fi <- df_ccf_iso %>% 
  filter(Alignment != "mis") %>% 
  select(-Alignment) %>% 
  set_sdif("Cue") %>% 
  set_sdif("Congruency", 3) %>% 
  set_sdif("SD") %>% 
  mutate(ConInc = if_else(Congruency=="con", -2/3,
                         if_else(Congruency %in% c("inc", "iso"), 1/3, NaN)),
         IncIso = if_else(Congruency=="iso", 2/3,
                         if_else(Congruency %in% c("inc", "con"), -1/3, NaN)),
         SD_C = if_else(SD=="same", -.5, if_else(SD=="diff", .5, NaN)),
         ConInc_SD = ConInc * SD_C,
         IncIso_SD = IncIso * SD_C)

# saveRDS(df_fi, file=here("jubail", "df_fi.rds"))

str(df_fi)
## tibble [107,135 × 14] (S3: tbl_df/tbl/data.frame)
##  $ Subject   : Factor w/ 480 levels "subj001","subj002",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ Cue       : Factor w/ 2 levels "top","bot": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "top" "bot"
##   .. .. ..$ : chr "2-1"
##  $ Congruency: Factor w/ 3 levels "con","inc","iso": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:3, 1:2] -0.667 0.333 0.333 -0.333 -0.333 ...
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:3] "con" "inc" "iso"
##   .. .. ..$ : chr [1:2] "2-1" "3-2"
##  $ SD        : Factor w/ 2 levels "same","diff": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "same" "diff"
##   .. .. ..$ : chr "2-1"
##  $ Correct   : logi [1:107135] TRUE TRUE TRUE TRUE TRUE TRUE ...
##  $ isSame    : num [1:107135] 1 1 1 1 1 1 1 1 1 0 ...
##  $ RT        : num [1:107135] 651 1399 986 819 1362 ...
##  $ Trial_num : num [1:107135] 56 70 76 90 113 114 122 146 147 151 ...
##  $ StimGroup : chr [1:107135] "F1" "M5" "F2" "M4" ...
##  $ ConInc    : num [1:107135] -0.667 -0.667 -0.667 -0.667 -0.667 ...
##  $ IncIso    : num [1:107135] -0.333 -0.333 -0.333 -0.333 -0.333 ...
##  $ SD_C      : num [1:107135] -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 ...
##  $ ConInc_SD : num [1:107135] 0.333 0.333 0.333 0.333 0.333 ...
##  $ IncIso_SD : num [1:107135] 0.167 0.167 0.167 0.167 0.167 ...

2 Holistic processing effects

Linear mixed-effects models are used to examine the holistic processing effects in the part-whole, standard composite face task, and complete composite face task.

2.1 Steps to obtain the optimal model

  1. If the maximal model did not converge, correlations between random effects were removed, making the zero-correlation-parameter (ZCP) model.
  2. Principal component analysis implemented with rePCA() function was then used to identify random effects that explained less than 0.1% of the total variances; they were removed from the ZCP model to make the reduced model.
  3. The extended model was built by adding back the correlations between random effects in the reduced model.
  4. If the extended model did not converge, the random effects that explained less than 1% of total variances were identified by rePCA() and removed to make the updated extended model; this step was iterated until an extended model converged.
  5. The converged extended model was then compared to the reduced model via anova() function and the model that explained the data better (with smaller Bayesian Information Criterion) was used as the optimal model.
  6. All follow-up analyses were performed on the optimal model.

2.2 Part-whole task

2.2.1 Accuracy

2.2.1.1 Maximal model

# thisfile <- here("jubail", "lmm_pw_acc_max.rds")
# 
# if (file.exists(thisfile)) {
#   lmm_pw_acc_max <- readRDS(thisfile)
# 
# } else {
#   
#   lmm_pw_acc_max <-
#     glmer(Correct ~ PW * Feature +
#             (PW * Feature | Subject) +
#             (PW * Feature | StimGroup),
#           df_pw,
#           family = binomial(link = "logit"),
#           control = glmerControl(optCtrl = list(maxfun = 1e7)))
#   
#   saveRDS(lmm_pw_acc_max, thisfile)
# }
# 
# summary(lmm_pw_acc_max)

2.2.1.2 Zero-correlation-parameter model

thisfile <- here("jubail", "lmm_pw_acc_zcp.rds")

if (file.exists(thisfile)) {
  lmm_pw_acc_zcp <- readRDS(thisfile)
  
} else {
  
  lmm_pw_acc_zcp <- 
    glmer(Correct ~ PW * Feature + 
            (PW_C + Feature_C + PW_Feature || Subject) +
            (PW_C + Feature_C + PW_Feature || StimGroup),
          df_pw,
          family = binomial(link = "logit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_pw_acc_zcp, thisfile)
}

summary(lmm_pw_acc_zcp)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( logit )
## Formula: Correct ~ PW * Feature + (PW_C + Feature_C + PW_Feature || Subject) +      (PW_C + Feature_C + PW_Feature || StimGroup)
##    Data: df_pw
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
##  35277.0  35381.0 -17626.5  35253.0    42868 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -6.5497  0.1945  0.3143  0.4625  2.0098 
## 
## Random effects:
##  Groups      Name        Variance  Std.Dev.
##  Subject     (Intercept) 4.239e-01 0.651106
##  Subject.1   PW_C        8.547e-02 0.292345
##  Subject.2   Feature_C   4.018e-01 0.633878
##  Subject.3   PW_Feature  3.699e-06 0.001923
##  StimGroup   (Intercept) 1.319e-01 0.363119
##  StimGroup.1 PW_C        4.323e-02 0.207918
##  StimGroup.2 Feature_C   4.858e-01 0.697006
##  StimGroup.3 PW_Feature  3.323e-01 0.576498
## Number of obs: 42880, groups:  Subject, 455; StimGroup, 12
## 
## Fixed effects:
##                  Estimate Std. Error z value Pr(>|z|)    
## (Intercept)       1.89707    0.11052  17.166  < 2e-16 ***
## PW2-1            -0.77077    0.06882 -11.200  < 2e-16 ***
## Feature2-1        0.93358    0.20593   4.533  5.8e-06 ***
## PW2-1:Feature2-1  0.19536    0.17704   1.103     0.27    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) PW2-1  Ftr2-1
## PW2-1       -0.018              
## Feature2-1   0.008 -0.005       
## PW2-1:Ft2-1 -0.003  0.039 -0.012
## optimizer (Nelder_Mead) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.00822005 (tol = 0.002, component 1)
summary(rePCA(lmm_pw_acc_zcp))
## $Subject
## Importance of components:
##                          [,1]   [,2]    [,3]     [,4]
## Standard deviation     0.6511 0.6339 0.29234 0.001923
## Proportion of Variance 0.4652 0.4410 0.09379 0.000000
## Cumulative Proportion  0.4652 0.9062 1.00000 1.000000
## 
## $StimGroup
## Importance of components:
##                          [,1]   [,2]   [,3]    [,4]
## Standard deviation     0.6970 0.5765 0.3631 0.20792
## Proportion of Variance 0.4891 0.3346 0.1328 0.04352
## Cumulative Proportion  0.4891 0.8237 0.9565 1.00000

(PW_Feature | Subject) will be removed.

2.2.1.3 Reduced model

thisfile <- here("jubail", "lmm_pw_acc_rdc.rds")

if (file.exists(thisfile)) {
  lmm_pw_acc_rdc <- readRDS(thisfile)
  
} else {
  
  lmm_pw_acc_rdc <- 
    glmer(Correct ~ PW * Feature + 
            (PW_C + Feature_C || Subject) + # + PW_Feature
            (PW_C + Feature_C + PW_Feature || StimGroup),
          df_pw,
          family = binomial(link = "logit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  saveRDS(lmm_pw_acc_rdc, thisfile)
}

summary(lmm_pw_acc_rdc)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( logit )
## Formula: Correct ~ PW * Feature + (PW_C + Feature_C || Subject) + (PW_C +      Feature_C + PW_Feature || StimGroup)
##    Data: df_pw
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
##  35275.0  35370.3 -17626.5  35253.0    42869 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -6.5498  0.1945  0.3143  0.4625  2.0098 
## 
## Random effects:
##  Groups      Name        Variance Std.Dev.
##  Subject     (Intercept) 0.42394  0.6511  
##  Subject.1   PW_C        0.08541  0.2923  
##  Subject.2   Feature_C   0.40188  0.6339  
##  StimGroup   (Intercept) 0.13207  0.3634  
##  StimGroup.1 PW_C        0.04321  0.2079  
##  StimGroup.2 Feature_C   0.48555  0.6968  
##  StimGroup.3 PW_Feature  0.33222  0.5764  
## Number of obs: 42880, groups:  Subject, 455; StimGroup, 12
## 
## Fixed effects:
##                  Estimate Std. Error z value Pr(>|z|)    
## (Intercept)       1.89688    0.11060  17.151  < 2e-16 ***
## PW2-1            -0.77086    0.06881 -11.203  < 2e-16 ***
## Feature2-1        0.93329    0.20586   4.534  5.8e-06 ***
## PW2-1:Feature2-1  0.19547    0.17681   1.106    0.269    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) PW2-1  Ftr2-1
## PW2-1       -0.018              
## Feature2-1   0.008 -0.005       
## PW2-1:Ft2-1 -0.004  0.039 -0.012

2.2.1.4 Extended model

thisfile <- here("jubail", "lmm_pw_acc_etd.rds")

if (file.exists(thisfile)) {
  lmm_pw_acc_etd <- readRDS(thisfile)
  
} else {
  # lmm_pw_acc_etd <- 
  #   glmer(Correct ~ PW * Feature + 
  #           (PW_C + Feature_C | Subject) + # + PW_Feature
  #           (PW_C + Feature_C + PW_Feature | StimGroup),
  #         df_pw,
  #         family = binomial(link = "logit"),
  #         control = glmerControl(optCtrl = list(maxfun = 1e7)))
  # the |grad| is about 0.0110027, so try other optimizer
  
  # use optimx instead
  lmm_pw_acc_etd <- 
    glmer(Correct ~ PW * Feature + 
            (PW_C + Feature_C | Subject) + # + PW_Feature
            (PW_C + Feature_C + PW_Feature | StimGroup),
          df_pw,
          family = binomial(link = "logit"),
          control = glmerControl(optimizer = "optimx", 
                                 optCtrl = list(method = "nlminb", 
                                                starttests = FALSE, 
                                                kkt = FALSE)))
  
  saveRDS(lmm_pw_acc_etd, thisfile)
}

summary(lmm_pw_acc_etd)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( logit )
## Formula: Correct ~ PW * Feature + (PW_C + Feature_C | Subject) + (PW_C +      Feature_C + PW_Feature | StimGroup)
##    Data: df_pw
## Control: glmerControl(optimizer = "optimx", optCtrl = list(method = "nlminb",      starttests = FALSE, kkt = FALSE))
## 
##      AIC      BIC   logLik deviance df.resid 
##  35166.0  35339.3 -17563.0  35126.0    42860 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -8.2217  0.1739  0.3128  0.4731  1.7376 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev. Corr             
##  Subject   (Intercept) 0.50670  0.7118                    
##            PW_C        0.11433  0.3381   -0.77            
##            Feature_C   0.47941  0.6924    0.57 -0.36      
##  StimGroup (Intercept) 0.13540  0.3680                    
##            PW_C        0.04638  0.2154   -0.48            
##            Feature_C   0.48417  0.6958   -0.12 -0.33      
##            PW_Feature  0.32037  0.5660    0.28 -0.32  0.03
## Number of obs: 42880, groups:  Subject, 455; StimGroup, 12
## 
## Fixed effects:
##                  Estimate Std. Error z value Pr(>|z|)    
## (Intercept)       1.94102    0.11279  17.209  < 2e-16 ***
## PW2-1            -0.87202    0.07253 -12.022  < 2e-16 ***
## Feature2-1        1.04188    0.20645   5.047  4.5e-07 ***
## PW2-1:Feature2-1  0.09829    0.17525   0.561    0.575    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) PW2-1  Ftr2-1
## PW2-1       -0.458              
## Feature2-1  -0.073 -0.299       
## PW2-1:Ft2-1  0.240 -0.206  0.009

2.2.1.5 Optimal model

anova(lmm_pw_acc_rdc, lmm_pw_acc_etd)

lmm_pw_acc_etd will be used as the optimal model.

lmm_pw_acc_opt <- lmm_pw_acc_etd

2.2.1.6 Effects of interest

Accuracy in each condition:

emm_pw_acc <- emmeans(lmm_pw_acc_opt, ~ PW)
summary(emm_pw_acc, type = "response")

Plot of the part-whole effect in accuracy:

emmip(lmm_pw_acc_opt, ~ PW, CIs = TRUE, type = "response") 

The part-whole effect in accuracy (differences between log odds ratio):

emm_pwe_acc <- contrast(emm_pw_acc, "pairwise")
# summary(emm_pwe_acc, side="=", level=.9, infer=TRUE)
summary(emm_pwe_acc, side=">", infer=TRUE)

The part-whole effect in accuracy (ratio of the odds):

summary(emm_pwe_acc, side=">", infer=TRUE, type="response") # exp(.872)
emmeans(lmm_pw_acc_opt, ~ PW + Feature, type = "response")
##  PW    Feature  prob      SE  df asymp.LCL asymp.UCL
##  whole mouth   0.868 0.02078 Inf     0.821     0.903
##  part  mouth   0.723 0.03177 Inf     0.657     0.781
##  whole eyes    0.947 0.00844 Inf     0.927     0.961
##  part  eyes    0.886 0.01454 Inf     0.854     0.912
## 
## Confidence level used: 0.95 
## Intervals are back-transformed from the logit scale

The part-whole effect for each feature:

contrast(emmeans(lmm_pw_acc_opt, ~ PW | Feature), "pairwise")
## Feature = mouth:
##  contrast     estimate    SE  df z.ratio p.value
##  whole - part    0.921 0.125 Inf   7.384  <.0001
## 
## Feature = eyes:
##  contrast     estimate    SE  df z.ratio p.value
##  whole - part    0.823 0.102 Inf   8.102  <.0001
## 
## Results are given on the log odds ratio (not the response) scale.

Plot for the part-whole effect for each feature:

emmip(lmm_pw_acc_opt, ~ PW | Feature, CIs = TRUE, type = "response") 

2.2.1.7 Plot

plot_pw_acc <- emm_pw_acc %>%
  summary(type="response") %>% 
  as_tibble() %>% 
  mutate(task = "PW") %>% 
  ggplot(aes(PW, prob, group=task)) +
  geom_point(size = 2, color = two_colors[1]) + # position = position_dodge(width = 0.1),
  geom_line(linetype="solid", color = two_colors[1], size = 0.8) +
  geom_errorbar(aes(ymin = asymp.LCL, ymax = asymp.UCL), size=1.5, width=0, 
                alpha = .6, # position = position_dodge(width = 0.1),
                color = two_colors[1],
                show.legend = F) + 
  coord_cartesian(ylim = c(.75, 1)) +  # set the limit for y axis c(0, 1100)
  labs(x = "Condition", y = "Accuracy") +  # set the names for main, x and y axises
  # geom_text(label = c("***", ""),
  #           color = "red",
  #           size = 6, nudge_y = 0.05, nudge_x = 0.5) + # add starts to the significant columns
  NULL
# ggsave(filename = "pw_acc.pdf", plot_pw_acc, width = 8, height = 6)
plot_pw_acc

2.2.2 Correct response times

2.2.2.1 Maximal model

# thisfile <- here("jubail", "lmm_pw_rt_max.rds")
# 
# if (file.exists(thisfile)) {
#   lmm_pw_rt_max <- readRDS(thisfile)
# 
# } else {
#   lmm_pw_rt_max <-
#     lmer(log(RT) ~ PW * Feature +
#            (PW * Feature | Subject) +
#            (PW * Feature | StimGroup),
#          filter(df_pw, Correct),
#          control = lmerControl(optCtrl = list(maxfun = 1e7)))
#   
#   saveRDS(lmm_pw_rt_max, thisfile)
# }
# 
# summary(lmm_pw_rt_max)

2.2.2.2 Zero-correlation-parameter model

thisfile <- here("jubail", "lmm_pw_rt_zcp.rds")

if (file.exists(thisfile)) {
  lmm_pw_rt_zcp <- readRDS(thisfile)
  
} else {
  
  lmm_pw_rt_zcp <- 
    lmer(log(RT) ~ PW * Feature + 
           (PW_C + Feature_C + PW_Feature || Subject) +
           (PW_C + Feature_C + PW_Feature || StimGroup),
         filter(df_pw, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_pw_rt_zcp, thisfile)
}

summary(lmm_pw_rt_zcp)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ PW * Feature + (PW_C + Feature_C + PW_Feature || Subject) +      (PW_C + Feature_C + PW_Feature || StimGroup)
##    Data: filter(df_pw, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 18371.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -6.0453 -0.6197 -0.1046  0.5048 11.3480 
## 
## Random effects:
##  Groups      Name        Variance Std.Dev.
##  Subject     (Intercept) 0.040248 0.20062 
##  Subject.1   PW_C        0.011982 0.10946 
##  Subject.2   Feature_C   0.010262 0.10130 
##  Subject.3   PW_Feature  0.018347 0.13545 
##  StimGroup   (Intercept) 0.002669 0.05166 
##  StimGroup.1 PW_C        0.001120 0.03347 
##  StimGroup.2 Feature_C   0.009391 0.09691 
##  StimGroup.3 PW_Feature  0.004479 0.06692 
##  Residual                0.089831 0.29972 
## Number of obs: 35403, groups:  Subject, 455; StimGroup, 12
## 
## Fixed effects:
##                  Estimate Std. Error       df t value Pr(>|t|)    
## (Intercept)       7.11075    0.01771 21.25714 401.578  < 2e-16 ***
## PW2-1             0.03203    0.01141 17.33540   2.807  0.01197 *  
## Feature2-1       -0.10017    0.02856 11.63890  -3.507  0.00452 ** 
## PW2-1:Feature2-1  0.04494    0.02134 13.21077   2.105  0.05494 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) PW2-1  Ftr2-1
## PW2-1        0.002              
## Feature2-1  -0.001 -0.001       
## PW2-1:Ft2-1 -0.001 -0.007  0.002
summary(rePCA(lmm_pw_rt_zcp))
## $Subject
## Importance of components:
##                          [,1]   [,2]   [,3]   [,4]
## Standard deviation     0.6694 0.4519 0.3652 0.3380
## Proportion of Variance 0.4979 0.2270 0.1482 0.1269
## Cumulative Proportion  0.4979 0.7248 0.8731 1.0000
## 
## $StimGroup
## Importance of components:
##                          [,1]   [,2]   [,3]    [,4]
## Standard deviation     0.3233 0.2233 0.1724 0.11167
## Proportion of Variance 0.5318 0.2536 0.1512 0.06344
## Cumulative Proportion  0.5318 0.7854 0.9366 1.00000

No random effects need to be removed. Therefore, we apply the extended model (i.e., the maximal model in this case) directly.

2.2.2.3 Extended model

thisfile <- here("jubail", "lmm_pw_rt_etd.rds")

if (file.exists(thisfile)) {
  lmm_pw_rt_etd <- readRDS(thisfile)
  
} else {
  lmm_pw_rt_etd <- 
    lmer(log(RT) ~ PW * Feature + 
           (PW_C + Feature_C + PW_Feature | Subject) +
           (PW_C + Feature_C + PW_Feature | StimGroup),
         filter(df_pw, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_pw_rt_etd, thisfile)
}

summary(lmm_pw_rt_etd)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ PW * Feature + (PW_C + Feature_C + PW_Feature | Subject) +      (PW_C + Feature_C + PW_Feature | StimGroup)
##    Data: filter(df_pw, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 18286.9
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -6.0929 -0.6206 -0.1062  0.5061 11.4120 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev. Corr             
##  Subject   (Intercept) 0.040164 0.20041                   
##            PW_C        0.012247 0.11066  -0.16            
##            Feature_C   0.010174 0.10087  -0.10  0.00      
##            PW_Feature  0.018837 0.13725   0.21 -0.56 -0.12
##  StimGroup (Intercept) 0.002655 0.05152                   
##            PW_C        0.001139 0.03375  -0.59            
##            Feature_C   0.009423 0.09707  -0.06  0.08      
##            PW_Feature  0.004529 0.06730   0.61 -0.68 -0.22
##  Residual              0.089828 0.29971                   
## Number of obs: 35403, groups:  Subject, 455; StimGroup, 12
## 
## Fixed effects:
##                  Estimate Std. Error       df t value Pr(>|t|)    
## (Intercept)       7.11073    0.01767 21.30825 402.473  < 2e-16 ***
## PW2-1             0.03167    0.01151 17.37547   2.752  0.01341 *  
## Feature2-1       -0.10035    0.02860 11.62636  -3.508  0.00452 ** 
## PW2-1:Feature2-1  0.04573    0.02147 13.24818   2.130  0.05243 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) PW2-1  Ftr2-1
## PW2-1       -0.458              
## Feature2-1  -0.057  0.066       
## PW2-1:Ft2-1  0.499 -0.603 -0.197
summary(rePCA(lmm_pw_rt_etd))
## $Subject
## Importance of components:
##                          [,1]   [,2]   [,3]    [,4]
## Standard deviation     0.6916 0.4957 0.3369 0.26245
## Proportion of Variance 0.5277 0.2711 0.1252 0.07599
## Cumulative Proportion  0.5277 0.7988 0.9240 1.00000
## 
## $StimGroup
## Importance of components:
##                          [,1]   [,2]    [,3]    [,4]
## Standard deviation     0.3336 0.2588 0.11738 0.07460
## Proportion of Variance 0.5632 0.3389 0.06974 0.02817
## Cumulative Proportion  0.5632 0.9021 0.97183 1.00000

2.2.2.4 Optimal model

anova(lmm_pw_rt_zcp, lmm_pw_rt_etd, refit=FALSE)

Use lmm_pw_rt_etd as the optimal model.

lmm_pw_rt_opt <- lmm_pw_rt_zcp

2.2.2.5 Effects of interest

Correct RT in each condition:

emm_pw_rt <- emmeans(lmm_pw_rt_opt, ~ PW)
summary(emm_pw_rt, type="response")

Plot for the each condition (RT):

emmip(lmm_pw_rt_opt, ~ PW, CIs = TRUE, type = "response") 

The part-whole effect in RT (differences between log RT):

emm_pwe_rt <- contrast(emm_pw_rt, "pairwise")
summary(emm_pwe_rt, side="<", infer=TRUE)

The part-whole effect in RT (ratio of the RT):

summary(emm_pwe_rt, side="<", infer=TRUE, type = "response") # exp(-0.032)

Correct RT in each condition including features:

emmeans(lmm_pw_rt_opt, ~ PW + Feature, type="response")
##  PW    Feature response   SE  df asymp.LCL asymp.UCL
##  whole mouth       1282 30.8 Inf      1223      1344
##  part  mouth       1294 31.2 Inf      1234      1357
##  whole eyes        1134 27.2 Inf      1082      1189
##  part  eyes        1197 28.8 Inf      1142      1255
## 
## Degrees-of-freedom method: asymptotic 
## Confidence level used: 0.95 
## Intervals are back-transformed from the log scale

The part-whole effect for each feature:

contrast(emmeans(lmm_pw_rt_opt, ~ PW | Feature), "pairwise")
## Feature = mouth:
##  contrast     estimate     SE  df z.ratio p.value
##  whole - part -0.00956 0.0157 Inf  -0.610  0.5421
## 
## Feature = eyes:
##  contrast     estimate     SE  df z.ratio p.value
##  whole - part -0.05450 0.0156 Inf  -3.500  0.0005
## 
## Degrees-of-freedom method: asymptotic 
## Results are given on the log (not the response) scale.

Plot for the conditions of each feature:

emmip(lmm_pw_rt_opt, ~ PW | Feature, CIs = TRUE, type = "response") 

2.2.2.6 Plot

plot_pw_rt <- emm_pw_rt %>% 
  summary(type = "response") %>% 
  as_tibble() %>% 
  mutate(task = "PW") %>% 
  ggplot(aes(y = response, x = PW, group = task)) +
  geom_point(size = 2, color = two_colors[1]) + # position = position_dodge(width = 0.1), 
  geom_line(linetype = "solid", # position = position_dodge(width = 0.1),
            color = two_colors[1],
            size = 0.8) +
  geom_errorbar(aes(ymin = asymp.LCL, ymax = asymp.UCL), size=1.5, width=0, 
                alpha = .6, # position = position_dodge(width = 0.1),
                color = two_colors[1],
                show.legend = F) + 
  coord_cartesian(ylim = c(1100, 1300)) +  # set the limit for y axis c(0, 1100)
  labs(x = "Condition", y = "Correct response times (ms)") +  # set the names for main, x and y axises
  # geom_text(label = c("***", ""),
  #           color = "red",
  #           size = 6, nudge_y = 75, nudge_x = 0.5) + # add starts to the significant columns
  NULL
# ggsave(filename = "pw_rt.pdf", plot_pw_rt, width = 8, height = 6)
plot_pw_rt

2.3 Standard composite face task

2.3.1 Accuracy

2.3.1.1 Maximal model

# thisfile <- here("jubail", "lmm_scf_acc_max.rds")
# 
# if (file.exists(thisfile)) {
#   lmm_scf_acc_max <- readRDS(thisfile)
# 
# } else {
#   lmm_scf_acc_max <-
#     glmer(Correct ~ Alignment +
#             (Alignment | Subject) +
#             (Alignment | StimGroup),
#           df_scf,
#           family = binomial(link = "logit"),
#           control = glmerControl(optCtrl = list(maxfun = 1e7)))
#   
#   saveRDS(lmm_scf_acc_max, thisfile)
# }
# 
# summary(lmm_scf_acc_max)

2.3.1.2 Zero-correlation-parameter model

thisfile <- here("jubail", "lmm_scf_acc_zcp.rds")

if (file.exists(thisfile)) {
  lmm_scf_acc_zcp <- readRDS(thisfile)
  
} else {
  
  lmm_scf_acc_zcp <- 
    glmer(Correct ~ Alignment + 
            (Ali_C || Subject) +
            (Ali_C || StimGroup),
          df_scf,
          family = binomial(link = "logit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_scf_acc_zcp, thisfile)
}

summary(lmm_scf_acc_zcp)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( logit )
## Formula: Correct ~ Alignment + (Ali_C || Subject) + (Ali_C || StimGroup)
##    Data: df_scf
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
##  21680.0  21730.8 -10834.0  21668.0    35642 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -7.9618  0.1743  0.2469  0.3526  4.2437 
## 
## Random effects:
##  Groups      Name        Variance Std.Dev.
##  Subject     (Intercept) 0.98598  0.9930  
##  Subject.1   Ali_C       0.24823  0.4982  
##  StimGroup   (Intercept) 0.04554  0.2134  
##  StimGroup.1 Ali_C       0.05901  0.2429  
## Number of obs: 35648, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   2.53689    0.08543  29.696   <2e-16 ***
## Alignment2-1  0.89375    0.09002   9.928   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr)
## Alignmnt2-1 0.040
summary(rePCA(lmm_scf_acc_zcp))
## $Subject
## Importance of components:
##                          [,1]   [,2]
## Standard deviation     0.9930 0.4982
## Proportion of Variance 0.7989 0.2011
## Cumulative Proportion  0.7989 1.0000
## 
## $StimGroup
## Importance of components:
##                          [,1]   [,2]
## Standard deviation     0.2429 0.2134
## Proportion of Variance 0.5645 0.4355
## Cumulative Proportion  0.5645 1.0000

No random effects need to be removed. Therefore, we apply the extended model (i.e., the maximal model in this case) directly.

2.3.1.3 Extended model

thisfile <- here("jubail", "lmm_scf_acc_etd.rds")

if (file.exists(thisfile)) {
  lmm_scf_acc_etd <- readRDS(thisfile)
  
} else {
  
  lmm_scf_acc_etd <- 
    glmer(Correct ~ Alignment + 
            (Ali_C | Subject) +
            (Ali_C | StimGroup),
          df_scf,
          family = binomial(link = "logit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  saveRDS(lmm_scf_acc_etd, thisfile)
}

summary(lmm_scf_acc_etd)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( logit )
## Formula: Correct ~ Alignment + (Ali_C | Subject) + (Ali_C | StimGroup)
##    Data: df_scf
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
##  21641.9  21709.7 -10812.9  21625.9    35640 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -8.7393  0.1608  0.2465  0.3578  3.5831 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev. Corr 
##  Subject   (Intercept) 1.05342  1.0264        
##            Ali_C       0.31534  0.5616   0.51 
##  StimGroup (Intercept) 0.04479  0.2116        
##            Ali_C       0.05357  0.2314   -1.00
## Number of obs: 35648, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   2.56810    0.08633   29.75   <2e-16 ***
## Alignment2-1  1.03098    0.09235   11.16   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr)
## Alignmnt2-1 -0.461
## optimizer (Nelder_Mead) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')
summary(rePCA(lmm_scf_acc_etd))
## $Subject
## Importance of components:
##                          [,1]   [,2]
## Standard deviation     1.0746 0.4625
## Proportion of Variance 0.8437 0.1563
## Cumulative Proportion  0.8437 1.0000
## 
## $StimGroup
## Importance of components:
##                          [,1] [,2]
## Standard deviation     0.3136    0
## Proportion of Variance 1.0000    0
## Cumulative Proportion  1.0000    1

By-item (Intercept) will be removed.

thisfile <- here("jubail", "lmm_scf_acc_etd2.rds")

if (file.exists(thisfile)) {
  
  lmm_scf_acc_etd2 <- readRDS(thisfile)
  
} else {
  
  lmm_scf_acc_etd2 <- 
    glmer(Correct ~ Alignment + 
            (Ali_C | Subject) +
            (0 + Ali_C | StimGroup),
          df_scf,
          family = binomial(link = "logit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  saveRDS(lmm_scf_acc_etd2, thisfile)
}

summary(lmm_scf_acc_etd2)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( logit )
## Formula: Correct ~ Alignment + (Ali_C | Subject) + (0 + Ali_C | StimGroup)
##    Data: df_scf
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
##  21750.5  21801.4 -10869.3  21738.5    35642 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -9.4869  0.1626  0.2530  0.3605  3.4895 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev. Corr
##  Subject   (Intercept) 1.0448   1.0222       
##            Ali_C       0.3208   0.5664   0.53
##  StimGroup Ali_C       0.1268   0.3560       
## Number of obs: 35648, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept)    2.5606     0.0543  47.156   <2e-16 ***
## Alignment2-1   1.0697     0.1257   8.509   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr)
## Alignmnt2-1 0.184

2.3.1.4 Optimal model

anova(lmm_scf_acc_zcp, lmm_scf_acc_etd2, refit=FALSE)

lmm_scf_acc_zcp will be used as the optimal model.

lmm_scf_acc_opt <- lmm_scf_acc_zcp

2.3.1.5 Effects of interest

Accuracy in each condition:

emm_scf_acc <- emmeans(lmm_scf_acc_opt, ~ Alignment)
summary(emm_scf_acc, type="response")

Plot of accuracy:

emmip(lmm_scf_acc_opt, ~ Alignment, CIs = TRUE, type = "response") 

The composite effect in accuracy (differences between log odds):

emm_scfe_acc <- contrast(emm_scf_acc, "pairwise")
summary(emm_scfe_acc, side="<", infer=TRUE)

The composite effect in accuracy (differences between log odds):

summary(emm_scfe_acc, side="<", infer=TRUE, type="response")

2.3.1.6 Plot

plot_scf_acc <- emm_scf_acc %>%
  summary(type="response") %>% 
  as_tibble() %>% 
  mutate(task = "SCF",
         Alignment = fct_recode(Alignment, aligned="ali", misaligned="mis")) %>% 
  ggplot(aes(Alignment, prob, group=task)) +
  geom_point(size = 2, color = two_colors[2]) + # position = position_dodge(width = 0.1),
  geom_line(linetype="dashed", color = two_colors[2], size = 0.8) +
  geom_errorbar(aes(ymin = asymp.LCL, ymax = asymp.UCL), size=1.5, width=0, 
                alpha = .6, # position = position_dodge(width = 0.1),
                color = two_colors[2],
                show.legend = F) + 
  coord_cartesian(ylim = c(.75, 1)) +  # set the limit for y axis c(0, 1100)
  labs(x = "Alignment", y = "Accuracy") +  # set the names for main, x and y axises
  # geom_text(label = c("***", ""),
  #           color = "red",
  #           size = 6, nudge_y = .08, nudge_x = 0.5) + # add starts to the significant columns
  NULL
# ggsave(filename = "scf_acc.pdf", plot_scf_acc, width = 8, height = 6)
plot_scf_acc

2.3.2 Correct response times

2.3.2.1 Maximal model

# thisfile <- here("jubail", "lmm_scf_rt_max.rds")
# 
# if (file.exists(thisfile)) {
#   lmm_scf_rt_max <- readRDS(thisfile)
# 
# } else {
#   
#   lmm_scf_rt_max <-
#     lmer(log(RT) ~ Alignment +
#            (Alignment | Subject) +
#            (Alignment | StimGroup),
#          filter(df_scf, Correct),
#          control = lmerControl(optCtrl = list(maxfun = 1e7)))
#   
#   saveRDS(lmm_scf_rt_max, thisfile)
# }
# 
# summary(lmm_scf_rt_max)

2.3.2.2 Zero-correlation-parameter model

thisfile <- here("jubail", "lmm_scf_rt_zcp.rds")

if (file.exists(thisfile)) {
  lmm_scf_rt_zcp <- readRDS(thisfile)
  
} else {
  
  lmm_scf_rt_zcp <- 
    lmer(log(RT) ~ Alignment + 
           (Ali_C || Subject) +
           (Ali_C || StimGroup),
         filter(df_scf, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_scf_rt_zcp, thisfile)
}

summary(lmm_scf_rt_zcp)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ Alignment + (Ali_C || Subject) + (Ali_C || StimGroup)
##    Data: filter(df_scf, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 8023.4
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -5.0501 -0.6423 -0.1114  0.5142 11.9073 
## 
## Random effects:
##  Groups      Name        Variance  Std.Dev.
##  Subject     (Intercept) 0.0410908 0.20271 
##  Subject.1   Ali_C       0.0022934 0.04789 
##  StimGroup   (Intercept) 0.0005342 0.02311 
##  StimGroup.1 Ali_C       0.0008913 0.02985 
##  Residual                0.0708967 0.26626 
## Number of obs: 31715, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##              Estimate Std. Error       df t value Pr(>|t|)    
## (Intercept)   6.52194    0.01209 58.70444 539.576  < 2e-16 ***
## Alignment2-1 -0.08054    0.01016  9.94442  -7.926 1.32e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr)
## Alignmnt2-1 -0.002
summary(rePCA(lmm_scf_rt_zcp))
## $Subject
## Importance of components:
##                          [,1]    [,2]
## Standard deviation     0.7613 0.17986
## Proportion of Variance 0.9471 0.05286
## Cumulative Proportion  0.9471 1.00000
## 
## $StimGroup
## Importance of components:
##                          [,1]   [,2]
## Standard deviation     0.1121 0.0868
## Proportion of Variance 0.6253 0.3747
## Cumulative Proportion  0.6253 1.0000

No random effects need to be removed. Therefore, we apply the extended model (i.e., the maximal model in this case) directly.

2.3.2.3 Extended model

thisfile <- here("jubail", "lmm_scf_rt_etd.rds")

if (file.exists(thisfile)) {
  lmm_scf_rt_etd <- readRDS(thisfile)
  
} else {
  
  lmm_scf_rt_etd <-
    lmer(log(RT) ~ Alignment +
           (Ali_C | Subject) +
           (Ali_C | StimGroup),
         filter(df_scf, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  saveRDS(lmm_scf_rt_etd, thisfile)
}

summary(lmm_scf_rt_etd)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ Alignment + (Ali_C | Subject) + (Ali_C | StimGroup)
##    Data: filter(df_scf, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 8011.3
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -5.0289 -0.6412 -0.1138  0.5145 11.9047 
## 
## Random effects:
##  Groups    Name        Variance  Std.Dev. Corr 
##  Subject   (Intercept) 0.0411141 0.20277       
##            Ali_C       0.0022733 0.04768  -0.17
##  StimGroup (Intercept) 0.0005369 0.02317       
##            Ali_C       0.0008972 0.02995  -0.81
##  Residual              0.0709008 0.26627       
## Number of obs: 31715, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##              Estimate Std. Error       df t value Pr(>|t|)    
## (Intercept)   6.52193    0.01210 58.41448 538.973  < 2e-16 ***
## Alignment2-1 -0.08082    0.01019  9.92938  -7.934 1.32e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr)
## Alignmnt2-1 -0.489

2.3.2.4 Optimal model

anova(lmm_scf_rt_zcp, lmm_scf_rt_etd, refit=FALSE)

lmm_scf_rt_etd will be used as the optimal model.

lmm_scf_rt_opt <- lmm_scf_rt_zcp

2.3.2.5 Effects of interest

Correct RT in each condition:

emm_scf_rt <- emmeans(lmm_scf_rt_opt, ~ Alignment)
summary(emm_scf_rt, type="response")

Plot of RT:

emmip(lmm_scf_rt_opt, ~ Alignment, CIs = TRUE, type = "response") 

The composite effect in RT (differences between log RT):

emm_scfe_rt <- contrast(emm_scf_rt, "pairwise")
summary(emm_scfe_rt, side=">", infer=TRUE)

The composite effect in RT (ratio of the RT):

summary(emm_scfe_rt, side=">", infer=TRUE, type = "response") # exp(0.0805)

2.3.2.6 Plot

plot_scf_rt <- emm_scf_rt %>% 
  summary(type = "response") %>% 
  as_tibble() %>% 
  mutate(task = "SCF",
         Alignment = fct_recode(Alignment, aligned="ali", misaligned="mis")) %>% 
  ggplot(aes(y = response, x = Alignment, group = task)) +
  geom_point(size = 2, color = two_colors[2]) + # position = position_dodge(width = 0.1), 
  geom_line(linetype = "dashed", # position = position_dodge(width = 0.1),
            color = two_colors[2],
            size = 0.8) +
  geom_errorbar(aes(ymin = asymp.LCL, ymax = asymp.UCL), size=1.5, width=0, 
                alpha = .6, # position = position_dodge(width = 0.1),
                color = two_colors[2],
                show.legend = F) + 
  coord_cartesian(ylim = c(600, 750)) +  # set the limit for y axis c(0, 1100)
  labs(x = "Alignment", y = "Correct response times (ms)") +  # set the names for main, x and y axises
  # geom_text(label = c("***", ""),
  #           color = "red",
  #           size = 6, nudge_y = 75, nudge_x = 0.5) + # add starts to the significant columns
  NULL
# ggsave(filename = "scf_rt.pdf", plot_scf_rt, width = 8, height = 6)
plot_scf_rt

2.4 Complete composite face task

2.4.1 Sensitivity d’

2.4.1.1 Maximal model

# thisfile <- here("jubail", "lmm_ccf_d_max.rds")
# 
# if (file.exists(thisfile)) {
#   lmm_ccf_d_max <- readRDS(thisfile)
#   
# } else {
#   
#   lmm_ccf_d_max <-
#     glmer(isSame ~ Congruency * Alignment * SD + Cue +
#             (Congruency * Alignment * SD | Subject) +
#             (Congruency * Alignment * SD | StimGroup),
#           df_ccf,
#           family = binomial(link = "probit"),
#           control = glmerControl(optCtrl = list(maxfun = 1e7)))
#   
#   saveRDS(lmm_ccf_d_max, thisfile)
# }
# 
# summary(lmm_ccf_d_max)

2.4.1.2 Zero-correlation-parameter model

thisfile <- here("jubail", "lmm_ccf_d_zcp.rds")

if (file.exists(thisfile)) {
  lmm_ccf_d_zcp <- readRDS(thisfile)
  
} else {
  
  lmm_ccf_d_zcp <- 
    glmer(isSame ~ Congruency * Alignment * SD + Cue +
            (Con_C + Ali_C + SD_C + 
               Con_Ali + Con_SD + Ali_SD + 
               Con_Ali_SD || Subject) +
            (Con_C + Ali_C + SD_C + 
               Con_Ali + Con_SD + Ali_SD + 
               Con_Ali_SD || StimGroup),
          df_ccf,
          family = binomial(link = "probit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_ccf_d_zcp, thisfile)
}

summary(lmm_ccf_d_zcp)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( probit )
## Formula: isSame ~ Congruency * Alignment * SD + Cue + (Con_C + Ali_C +      SD_C + Con_Ali + Con_SD + Ali_SD + Con_Ali_SD || Subject) +      (Con_C + Ali_C + SD_C + Con_Ali + Con_SD + Ali_SD + Con_Ali_SD ||          StimGroup)
##    Data: df_ccf
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
## 141572.2 141819.0 -70761.1 141522.2   142977 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -7.1001 -0.5842  0.2328  0.5105  5.3709 
## 
## Random effects:
##  Groups      Name        Variance  Std.Dev. 
##  Subject     (Intercept) 6.598e-02 0.2568624
##  Subject.1   Con_C       2.716e-04 0.0164797
##  Subject.2   Ali_C       3.111e-02 0.1763923
##  Subject.3   SD_C        2.908e-01 0.5392848
##  Subject.4   Con_Ali     3.693e-07 0.0006077
##  Subject.5   Con_SD      1.419e-01 0.3766749
##  Subject.6   Ali_SD      2.073e-04 0.0143983
##  Subject.7   Con_Ali_SD  6.819e-02 0.2611273
##  StimGroup   (Intercept) 9.431e-03 0.0971136
##  StimGroup.1 Con_C       3.631e-04 0.0190545
##  StimGroup.2 Ali_C       2.896e-04 0.0170190
##  StimGroup.3 SD_C        2.641e-02 0.1624970
##  StimGroup.4 Con_Ali     1.048e-03 0.0323716
##  StimGroup.5 Con_SD      3.468e-03 0.0588870
##  StimGroup.6 Ali_SD      7.814e-04 0.0279535
##  StimGroup.7 Con_Ali_SD  4.750e-07 0.0006892
## Number of obs: 143002, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                                  Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                       0.17361    0.03323   5.225 1.74e-07 ***
## Congruency2-1                    -0.10965    0.01001 -10.951  < 2e-16 ***
## Alignment2-1                      0.10676    0.01265   8.442  < 2e-16 ***
## SD2-1                            -1.65282    0.05782 -28.587  < 2e-16 ***
## Cue2-1                            0.12292    0.00773  15.903  < 2e-16 ***
## Congruency2-1:Alignment2-1        0.16714    0.01880   8.891  < 2e-16 ***
## Congruency2-1:SD2-1               0.78710    0.03020  26.065  < 2e-16 ***
## Alignment2-1:SD2-1               -0.04675    0.01815  -2.576     0.01 ** 
## Congruency2-1:Alignment2-1:SD2-1 -0.64145    0.03379 -18.982  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Alg2-1 SD2-1  Cue2-1 Cn2-1:A2-1 C2-1:S A2-1:S
## Congrncy2-1 -0.013                                                     
## Alignmnt2-1  0.000  0.025                                              
## SD2-1       -0.003  0.009 -0.002                                       
## Cue2-1       0.000  0.000  0.002 -0.005                                
## Cng2-1:A2-1  0.005 -0.004 -0.066 -0.004 -0.005                         
## Cn2-1:SD2-1  0.005 -0.053 -0.010 -0.011  0.005 -0.007                  
## Al2-1:SD2-1 -0.002 -0.021 -0.067  0.000  0.004  0.048      0.023       
## C2-1:A2-1:S -0.004 -0.012  0.038  0.007 -0.002 -0.093     -0.004 -0.102
## optimizer (Nelder_Mead) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.0242223 (tol = 0.002, component 1)
summary(rePCA(lmm_ccf_d_zcp))
## $Subject
## Importance of components:
##                          [,1]   [,2]   [,3]   [,4]    [,5]    [,6]    [,7]      [,8]
## Standard deviation     0.5393 0.3767 0.2611 0.2569 0.17639 0.01648 0.01440 0.0006077
## Proportion of Variance 0.4859 0.2371 0.1139 0.1102 0.05199 0.00045 0.00035 0.0000000
## Cumulative Proportion  0.4859 0.7230 0.8370 0.9472 0.99920 0.99965 1.00000 1.0000000
## 
## $StimGroup
## Importance of components:
##                          [,1]    [,2]    [,3]    [,4]    [,5]    [,6]    [,7]      [,8]
## Standard deviation     0.1625 0.09711 0.05889 0.03237 0.02795 0.01905 0.01702 0.0006892
## Proportion of Variance 0.6319 0.22570 0.08299 0.02508 0.01870 0.00869 0.00693 0.0000100
## Cumulative Proportion  0.6319 0.85760 0.94059 0.96567 0.98437 0.99306 0.99999 1.0000000

By-subject Con_Ali, Con_C, Ali_SD, and by-item Con_Ali_SD were removed.

2.4.1.3 Reduced model

thisfile <- here("jubail", "lmm_ccf_d_rdc.rds")

if (file.exists(thisfile)) {
  lmm_ccf_d_rdc <- readRDS(thisfile)
  
} else {
  
  lmm_ccf_d_rdc <- 
    glmer(isSame ~ Congruency * Alignment * SD + Cue +
            (Ali_C + SD_C + # Con_C + 
               Con_SD +  # Con_Ali + Ali_SD +
               Con_Ali_SD || Subject) +
            (Con_C + Ali_C + SD_C + 
               Con_Ali + Con_SD + Ali_SD # + 
             || StimGroup), # Con_Ali_SD
          df_ccf,
          family = binomial(link = "probit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_ccf_d_rdc, thisfile)
}

summary(lmm_ccf_d_rdc)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( probit )
## Formula: isSame ~ Congruency * Alignment * SD + Cue + (Ali_C + SD_C +      Con_SD + Con_Ali_SD || Subject) + (Con_C + Ali_C + SD_C +      Con_Ali + Con_SD + Ali_SD || StimGroup)
##    Data: df_ccf
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
## 141564.3 141771.6 -70761.1 141522.3   142981 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -7.1035 -0.5841  0.2329  0.5105  5.3774 
## 
## Random effects:
##  Groups      Name        Variance  Std.Dev.
##  Subject     (Intercept) 0.0659269 0.25676 
##  Subject.1   Ali_C       0.0311195 0.17641 
##  Subject.2   SD_C        0.2907778 0.53924 
##  Subject.3   Con_SD      0.1419338 0.37674 
##  Subject.4   Con_Ali_SD  0.0683220 0.26138 
##  StimGroup   (Intercept) 0.0094380 0.09715 
##  StimGroup.1 Con_C       0.0003623 0.01904 
##  StimGroup.2 Ali_C       0.0002877 0.01696 
##  StimGroup.3 SD_C        0.0263888 0.16245 
##  StimGroup.4 Con_Ali     0.0010548 0.03248 
##  StimGroup.5 Con_SD      0.0034791 0.05898 
##  StimGroup.6 Ali_SD      0.0007764 0.02786 
## Number of obs: 143002, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                                   Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                       0.173816   0.033236   5.230  1.7e-07 ***
## Congruency2-1                    -0.109583   0.009964 -10.998  < 2e-16 ***
## Alignment2-1                      0.106833   0.012639   8.453  < 2e-16 ***
## SD2-1                            -1.652576   0.057787 -28.598  < 2e-16 ***
## Cue2-1                            0.122877   0.007729  15.898  < 2e-16 ***
## Congruency2-1:Alignment2-1        0.167076   0.018815   8.880  < 2e-16 ***
## Congruency2-1:SD2-1               0.786786   0.030206  26.048  < 2e-16 ***
## Alignment2-1:SD2-1               -0.046773   0.018114  -2.582  0.00982 ** 
## Congruency2-1:Alignment2-1:SD2-1 -0.641429   0.033804 -18.975  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Alg2-1 SD2-1  Cue2-1 Cn2-1:A2-1 C2-1:S A2-1:S
## Congrncy2-1 -0.013                                                     
## Alignmnt2-1  0.000  0.025                                              
## SD2-1       -0.002  0.008 -0.002                                       
## Cue2-1       0.000  0.001  0.002 -0.005                                
## Cng2-1:A2-1  0.005 -0.004 -0.066 -0.004 -0.005                         
## Cn2-1:SD2-1  0.005 -0.052 -0.010 -0.011  0.005 -0.007                  
## Al2-1:SD2-1 -0.002 -0.020 -0.067  0.000  0.004  0.048      0.023       
## C2-1:A2-1:S -0.004 -0.012  0.038  0.007 -0.002 -0.093     -0.004 -0.102
summary(rePCA(lmm_ccf_d_rdc))
## $Subject
## Importance of components:
##                          [,1]   [,2]   [,3]   [,4]    [,5]
## Standard deviation     0.5392 0.3767 0.2614 0.2568 0.17641
## Proportion of Variance 0.4862 0.2373 0.1142 0.1102 0.05203
## Cumulative Proportion  0.4862 0.7235 0.8377 0.9480 1.00000
## 
## $StimGroup
## Importance of components:
##                          [,1]    [,2]    [,3]    [,4]    [,5]    [,6]    [,7]
## Standard deviation     0.1624 0.09715 0.05898 0.03248 0.02786 0.01904 0.01696
## Proportion of Variance 0.6315 0.22586 0.08326 0.02524 0.01858 0.00867 0.00688
## Cumulative Proportion  0.6315 0.85737 0.94062 0.96587 0.98444 0.99312 1.00000

2.4.1.4 Extended model

thisfile <- here("jubail", "lmm_ccf_d_etd.rds")

if (file.exists(thisfile)) {
  lmm_ccf_d_etd <- readRDS(thisfile)
  
} else {
  lmm_ccf_d_etd <-
    glmer(isSame ~ Congruency * Alignment * SD + Cue +
            (Ali_C + SD_C + # Con_C + 
               Con_SD +  # Con_Ali + Ali_SD +
               Con_Ali_SD | Subject) +
            (Con_C + Ali_C + SD_C + 
               Con_Ali + Con_SD + Ali_SD # + 
             | StimGroup), # Con_Ali_SD
          df_ccf,
          family = binomial(link = "probit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_ccf_d_etd, thisfile)
}

summary(lmm_ccf_d_etd)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( probit )
## Formula: isSame ~ Congruency * Alignment * SD + Cue + (Ali_C + SD_C +      Con_SD + Con_Ali_SD | Subject) + (Con_C + Ali_C + SD_C +      Con_Ali + Con_SD + Ali_SD | StimGroup)
##    Data: df_ccf
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
## 141470.9 141984.2 -70683.5 141366.9   142950 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -7.4641 -0.5858  0.2267  0.5143  5.7349 
## 
## Random effects:
##  Groups    Name        Variance  Std.Dev. Corr                               
##  Subject   (Intercept) 0.0663148 0.25752                                     
##            Ali_C       0.0311955 0.17662   0.06                              
##            SD_C        0.2970462 0.54502  -0.25 -0.17                        
##            Con_SD      0.1431106 0.37830  -0.03  0.09 -0.31                  
##            Con_Ali_SD  0.0686201 0.26195  -0.09 -0.09  0.64 -0.93            
##  StimGroup (Intercept) 0.0089202 0.09445                                     
##            Con_C       0.0008627 0.02937   0.97                              
##            Ali_C       0.0004587 0.02142   0.15 -0.08                        
##            SD_C        0.0267054 0.16342   0.97  0.94  0.14                  
##            Con_Ali     0.0022944 0.04790  -0.90 -0.83 -0.33 -0.97            
##            Con_SD      0.0048947 0.06996  -0.77 -0.85  0.33 -0.63  0.44      
##            Ali_SD      0.0009042 0.03007   0.21  0.31 -0.43  0.36 -0.28 -0.04
## Number of obs: 143002, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                                  Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                       0.17605    0.03244   5.427 5.72e-08 ***
## Congruency2-1                    -0.11074    0.01228  -9.015  < 2e-16 ***
## Alignment2-1                      0.10732    0.01336   8.035 9.39e-16 ***
## SD2-1                            -1.65567    0.05815 -28.470  < 2e-16 ***
## Cue2-1                            0.12287    0.00773  15.894  < 2e-16 ***
## Congruency2-1:Alignment2-1        0.17021    0.02206   7.716 1.20e-14 ***
## Congruency2-1:SD2-1               0.79575    0.03256  24.441  < 2e-16 ***
## Alignment2-1:SD2-1               -0.03973    0.01867  -2.128   0.0333 *  
## Congruency2-1:Alignment2-1:SD2-1 -0.66219    0.03417 -19.378  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Alg2-1 SD2-1  Cue2-1 Cn2-1:A2-1 C2-1:S A2-1:S
## Congrncy2-1  0.664                                                     
## Alignmnt2-1  0.084 -0.008                                              
## SD2-1        0.744  0.637  0.018                                       
## Cue2-1       0.000  0.000  0.001 -0.005                                
## Cng2-1:A2-1 -0.559 -0.433 -0.177 -0.592 -0.004                         
## Cn2-1:SD2-1 -0.482 -0.481  0.134 -0.467  0.005  0.200                  
## Al2-1:SD2-1  0.095  0.100 -0.180  0.161  0.005 -0.055      0.013       
## C2-1:A2-1:S -0.017 -0.010  0.019  0.109 -0.002 -0.090     -0.187 -0.110
## optimizer (Nelder_Mead) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')
summary(rePCA(lmm_ccf_d_etd))
## $Subject
## Importance of components:
##                          [,1]   [,2]    [,3]    [,4]    [,5]
## Standard deviation     0.6104 0.3795 0.24355 0.17409 0.00681
## Proportion of Variance 0.6146 0.2375 0.09784 0.04999 0.00008
## Cumulative Proportion  0.6146 0.8521 0.94993 0.99992 1.00000
## 
## $StimGroup
## Importance of components:
##                          [,1]    [,2]    [,3]    [,4]      [,5]      [,6]      [,7]
## Standard deviation     0.2012 0.05752 0.03331 0.01247 0.0001238 7.184e-05 4.581e-06
## Proportion of Variance 0.8984 0.07347 0.02464 0.00346 0.0000000 0.000e+00 0.000e+00
## Cumulative Proportion  0.8984 0.97190 0.99654 1.00000 1.0000000 1.000e+00 1.000e+00

By-subject Ali_C, by-item Ali_C, Con_C, Ali_SD, and Con_Ali were removed.

thisfile <- here("jubail", "lmm_ccf_d_etd2.rds")

if (file.exists(thisfile)) {
  lmm_ccf_d_etd2 <- readRDS(thisfile)
  
} else {
  
  lmm_ccf_d_etd2 <-
    glmer(isSame ~ Congruency * Alignment * SD + Cue +
            (SD_C + # Con_C + Ali_C + 
               Con_SD +  # Con_Ali + Ali_SD +
               Con_Ali_SD | Subject) +
            (SD_C + # Con_C + Ali_C + 
               Con_SD  # + + Ali_SD Con_Ali + 
             | StimGroup), # Con_Ali_SD
          df_ccf,
          family = binomial(link = "probit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_ccf_d_etd2, thisfile)
}

summary(lmm_ccf_d_etd2)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( probit )
## Formula: isSame ~ Congruency * Alignment * SD + Cue + (SD_C + Con_SD +      Con_Ali_SD | Subject) + (SD_C + Con_SD | StimGroup)
##    Data: df_ccf
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
## 141625.1 141871.9 -70787.6 141575.1   142977 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -7.8670 -0.5889  0.2314  0.5162  6.1525 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev. Corr             
##  Subject   (Intercept) 0.065667 0.25626                   
##            SD_C        0.294254 0.54245  -0.25            
##            Con_SD      0.142578 0.37759  -0.03 -0.32      
##            Con_Ali_SD  0.070031 0.26463  -0.08  0.65 -0.93
##  StimGroup (Intercept) 0.009245 0.09615                   
##            SD_C        0.025594 0.15998   0.97            
##            Con_SD      0.003538 0.05948  -0.81 -0.63      
## Number of obs: 143002, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                                   Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                       0.175162   0.032930   5.319 1.04e-07 ***
## Congruency2-1                    -0.109822   0.007997 -13.733  < 2e-16 ***
## Alignment2-1                      0.103077   0.007845  13.139  < 2e-16 ***
## SD2-1                            -1.647772   0.057169 -28.823  < 2e-16 ***
## Cue2-1                            0.122396   0.007714  15.866  < 2e-16 ***
## Congruency2-1:Alignment2-1        0.170973   0.015895  10.757  < 2e-16 ***
## Congruency2-1:SD2-1               0.795063   0.030368  26.181  < 2e-16 ***
## Alignment2-1:SD2-1               -0.034321   0.015745  -2.180   0.0293 *  
## Congruency2-1:Alignment2-1:SD2-1 -0.667571   0.034094 -19.580  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Alg2-1 SD2-1  Cue2-1 Cn2-1:A2-1 C2-1:S A2-1:S
## Congrncy2-1 -0.018                                                     
## Alignmnt2-1 -0.001  0.055                                              
## SD2-1        0.747  0.011 -0.002                                       
## Cue2-1       0.000  0.001  0.002 -0.005                                
## Cng2-1:A2-1  0.008 -0.011 -0.130 -0.005 -0.005                         
## Cn2-1:SD2-1 -0.462 -0.070 -0.017 -0.441  0.005 -0.007                  
## Al2-1:SD2-1 -0.002 -0.031 -0.122 -0.002  0.006  0.067      0.030       
## C2-1:A2-1:S -0.015 -0.013  0.063  0.113 -0.002 -0.122     -0.203 -0.122
## optimizer (Nelder_Mead) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.0296896 (tol = 0.002, component 1)
summary(rePCA(lmm_ccf_d_etd2))
## $Subject
## Importance of components:
##                          [,1]   [,2]   [,3]     [,4]
## Standard deviation     0.6095 0.3774 0.2422 0.005865
## Proportion of Variance 0.6488 0.2487 0.1024 0.000060
## Cumulative Proportion  0.6488 0.8975 0.9999 1.000000
## 
## $StimGroup
## Importance of components:
##                          [,1]    [,2]     [,3]
## Standard deviation     0.1901 0.04726 0.000511
## Proportion of Variance 0.9418 0.05821 0.000010
## Cumulative Proportion  0.9418 0.99999 1.000000

By-subject (Intercept) and by-item Con_SD were removed.

thisfile <- here("jubail", "lmm_ccf_d_etd3.rds")

if (file.exists(thisfile)) {
  lmm_ccf_d_etd3 <- readRDS(thisfile)
  
} else {
  
  lmm_ccf_d_etd3 <-
    glmer(isSame ~ Congruency * Alignment * SD + Cue +
            (0 + SD_C + # Con_C + Ali_C + 
               Con_SD +  # Con_Ali + Ali_SD +
               Con_Ali_SD | Subject) +
            (SD_C # +  Con_C + Ali_C + 
             # + Con_SD + Ali_SD Con_Ali + 
             | StimGroup), # Con_Ali_SD
          df_ccf,
          family = binomial(link = "probit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_ccf_d_etd3, thisfile)
}

summary(lmm_ccf_d_etd3)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( probit )
## Formula: isSame ~ Congruency * Alignment * SD + Cue + (0 + SD_C + Con_SD +      Con_Ali_SD | Subject) + (SD_C | StimGroup)
##    Data: df_ccf
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
## 144820.3 144998.0 -72392.2 144784.3   142984 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -6.6749 -0.6002  0.2810  0.5243  4.7748 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev. Corr       
##  Subject   SD_C        0.27471  0.52413             
##            Con_SD      0.13864  0.37235  -0.33      
##            Con_Ali_SD  0.06937  0.26339   0.63 -0.94
##  StimGroup (Intercept) 0.00900  0.09487             
##            SD_C        0.02417  0.15548  0.97       
## Number of obs: 143002, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                                   Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                       0.160974   0.030235   5.324 1.02e-07 ***
## Congruency2-1                    -0.103889   0.007718 -13.461  < 2e-16 ***
## Alignment2-1                      0.099515   0.007709  12.910  < 2e-16 ***
## SD2-1                            -1.592282   0.055496 -28.692  < 2e-16 ***
## Cue2-1                            0.118102   0.007600  15.539  < 2e-16 ***
## Congruency2-1:Alignment2-1        0.161667   0.015425  10.481  < 2e-16 ***
## Congruency2-1:SD2-1               0.764995   0.023452  32.619  < 2e-16 ***
## Alignment2-1:SD2-1               -0.030491   0.015482  -1.969   0.0489 *  
## Congruency2-1:Alignment2-1:SD2-1 -0.641519   0.033546 -19.124  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Alg2-1 SD2-1  Cue2-1 Cn2-1:A2-1 C2-1:S A2-1:S
## Congrncy2-1 -0.016                                                     
## Alignmnt2-1 -0.001  0.052                                              
## SD2-1        0.846  0.009 -0.002                                       
## Cue2-1      -0.001  0.001  0.002 -0.005                                
## Cng2-1:A2-1  0.007 -0.008 -0.121 -0.004 -0.006                         
## Cn2-1:SD2-1  0.006 -0.075 -0.020 -0.121  0.006 -0.010                  
## Al2-1:SD2-1 -0.002 -0.029 -0.108 -0.002  0.006  0.063      0.038       
## C2-1:A2-1:S -0.004 -0.014  0.059  0.111 -0.002 -0.103     -0.265 -0.118
## optimizer (Nelder_Mead) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.00256836 (tol = 0.002, component 1)
summary(rePCA(lmm_ccf_d_etd3))
## $Subject
## Importance of components:
##                          [,1]  [,2]     [,3]
## Standard deviation     0.5912 0.365 0.000185
## Proportion of Variance 0.7240 0.276 0.000000
## Cumulative Proportion  0.7240 1.000 1.000000
## 
## $StimGroup
## Importance of components:
##                          [,1]    [,2]
## Standard deviation     0.1809 0.02112
## Proportion of Variance 0.9866 0.01345
## Cumulative Proportion  0.9866 1.00000

By-subject Con_Ali_SD was removed.

thisfile <- here("jubail", "lmm_ccf_d_etd4.rds")

if (file.exists(thisfile)) {
  lmm_ccf_d_etd4 <- readRDS(thisfile)
  
} else {
  lmm_ccf_d_etd4 <-
    glmer(isSame ~ Congruency * Alignment * SD + Cue +
            (0 + SD_C + # Con_C + Ali_C + 
               Con_SD  # Con_Ali + Ali_SD + + 
             | Subject) + # Con_Ali_SD
            (SD_C # +  Con_C + Ali_C + 
             # + Con_SD + Ali_SD Con_Ali + 
             | StimGroup), # Con_Ali_SD
          df_ccf,
          family = binomial(link = "probit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_ccf_d_etd4, thisfile)
}

summary(lmm_ccf_d_etd4)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( probit )
## Formula: isSame ~ Congruency * Alignment * SD + Cue + (0 + SD_C + Con_SD |      Subject) + (SD_C | StimGroup)
##    Data: df_ccf
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
## 144867.4 145015.5 -72418.7 144837.4   142987 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -6.2945 -0.5994  0.2826  0.5229  4.4249 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev. Corr 
##  Subject   SD_C        0.273059 0.52255       
##            Con_SD      0.137954 0.37142  -0.33
##  StimGroup (Intercept) 0.008998 0.09486       
##            SD_C        0.024125 0.15532  0.97 
## Number of obs: 143002, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                                   Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                       0.160739   0.030238   5.316 1.06e-07 ***
## Congruency2-1                    -0.104048   0.007714 -13.488  < 2e-16 ***
## Alignment2-1                      0.100920   0.007703  13.102  < 2e-16 ***
## SD2-1                            -1.590898   0.055433 -28.700  < 2e-16 ***
## Cue2-1                            0.117995   0.007598  15.529  < 2e-16 ***
## Congruency2-1:Alignment2-1        0.158500   0.015405  10.289  < 2e-16 ***
## Congruency2-1:SD2-1               0.764075   0.023415  32.632  < 2e-16 ***
## Alignment2-1:SD2-1               -0.040905   0.015406  -2.655  0.00793 ** 
## Congruency2-1:Alignment2-1:SD2-1 -0.615645   0.030813 -19.980  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Alg2-1 SD2-1  Cue2-1 Cn2-1:A2-1 C2-1:S A2-1:S
## Congrncy2-1 -0.016                                                     
## Alignmnt2-1 -0.001  0.048                                              
## SD2-1        0.847  0.009 -0.002                                       
## Cue2-1      -0.001  0.001  0.002 -0.005                                
## Cng2-1:A2-1  0.006 -0.005 -0.121 -0.004 -0.006                         
## Cn2-1:SD2-1  0.006 -0.075 -0.019 -0.120  0.006 -0.011                  
## Al2-1:SD2-1 -0.002 -0.028 -0.106 -0.001  0.006  0.061      0.033       
## C2-1:A2-1:S -0.004 -0.015  0.061  0.007 -0.002 -0.106     -0.004 -0.121
summary(rePCA(lmm_ccf_d_etd4))
## $Subject
## Importance of components:
##                          [,1]   [,2]
## Standard deviation     0.5461 0.3359
## Proportion of Variance 0.7255 0.2745
## Cumulative Proportion  0.7255 1.0000
## 
## $StimGroup
## Importance of components:
##                          [,1]    [,2]
## Standard deviation     0.1808 0.02109
## Proportion of Variance 0.9866 0.01343
## Cumulative Proportion  0.9866 1.00000

2.4.1.5 Optimal model

anova(lmm_ccf_d_rdc, lmm_ccf_d_etd4, refit=FALSE)

lmm_ccf_d_rdc is used as the optimal model.

lmm_ccf_d_opt <- lmm_ccf_d_rdc

2.4.1.6 Effects of interest

Hit and false alarm rates:

(emm_ccf_probit <- emmeans(lmm_ccf_d_opt, ~ Congruency + Alignment + SD, type = "response"))
##  Congruency Alignment SD    prob      SE  df asymp.LCL asymp.UCL
##  con        ali       same 0.905 0.00790 Inf     0.888     0.919
##  inc        ali       same 0.713 0.01561 Inf     0.681     0.743
##  con        mis       same 0.884 0.00907 Inf     0.865     0.901
##  inc        mis       same 0.825 0.01187 Inf     0.801     0.848
##  con        ali       diff 0.191 0.01253 Inf     0.167     0.216
##  inc        ali       diff 0.304 0.01601 Inf     0.273     0.336
##  con        mis       diff 0.238 0.01419 Inf     0.211     0.266
##  inc        mis       diff 0.306 0.01606 Inf     0.275     0.338
## 
## Results are averaged over the levels of: Cue 
## Confidence level used: 0.95 
## Intervals are back-transformed from the probit scale

Sensitivity d’ in each condition:

emm_ccf_d <- contrast(emm_ccf_probit, "pairwise", simple="SD")
summary(emm_ccf_d[1:4], infer=TRUE)

Plot of d’ in each condition:

emmip(emm_ccf_d, Congruency ~ Alignment, CIs = TRUE)

Congruency effect of d’ for the aligned condition:

emm_con_d <- contrast(emm_ccf_probit, interaction="pairwise", by="Alignment")
summary(emm_con_d[1], side=">", infer=TRUE)

Composite effect of d’:

emm_ccfe_d <- contrast(emm_ccf_probit, interaction="pairwise")
summary(emm_ccfe_d, side=">", infer=TRUE)
# facilitation and interference of d'
contrast(emm_ccf_probit, interaction="pairwise", by="Congruency")
## Congruency = con:
##  Alignment_pairwise SD_pairwise estimate     SE  df z.ratio p.value
##  ali - mis          same - diff    0.274 0.0260 Inf  10.533  <.0001
## 
## Congruency = inc:
##  Alignment_pairwise SD_pairwise estimate     SE  df z.ratio p.value
##  ali - mis          same - diff   -0.367 0.0235 Inf -15.652  <.0001
## 
## Results are averaged over the levels of: Cue

2.4.1.7 Plot

plot_ccf_d <- emm_ccf_d %>%
  summary(infer=TRUE) %>% 
  as_tibble() %>% 
  mutate(Congruency = fct_recode(Congruency, congruent="con", incongruent="inc"),
         Alignment = fct_recode(Alignment, aligned="ali", misaligned="mis")) %>%
  ggplot(aes(Alignment, estimate, color=Congruency, group=Congruency)) +
  geom_point(size = 2) + # position = position_dodge(width = 0.1),
  geom_line(aes(linetype = Congruency), size = 0.8) +
  scale_linetype_manual(values=c("solid", "dashed")) +
  scale_color_manual(values=two_colors) +
  geom_errorbar(aes(ymin = asymp.LCL, ymax = asymp.UCL), size=1.5, width=0, 
                alpha = .6, # position = position_dodge(width = 0.1),
                show.legend = F) + 
  coord_cartesian(ylim = c(0,3.2)) +  # set the limit for y axis c(0, 1100)
  labs(y = expression("Sensitivity"~italic("d'")), fill = "Congruency") +  # set the names for main, x and y axises
  # geom_text(label = c("***", "", "", ""), 
  #           color = "red", 
  #           size = 6, nudge_y = 0.5, nudge_x = 0.5) + # add starts to the significant columns
  NULL
# ggsave(filename = "ccf_d.pdf", plot_ccf_d, width = 8, height = 6)
plot_ccf_d

plot_ccf_d_fi <- contrast(emm_ccf_probit, interaction="pairwise", by="Congruency") %>% 
  summary(infer=TRUE) %>% 
  as_tibble() %>% 
  mutate(Congruency = fct_recode(Congruency, congruent="con", incongruent="inc")) %>% 
  ggplot(aes(y = estimate, x = Congruency, color = Congruency)) +
  geom_point(size = 2) +
  geom_errorbar(aes(ymin = asymp.LCL, ymax = asymp.UCL), size=1.5, width=0, 
                alpha = .6) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  scale_color_manual(values=two_colors) +
  # coord_cartesian(ylim = ylimit_cf_fi_d) +  # set the limit for y axis c(0, 1100)
  labs(x = "Congruency", y = expression(italic("d'")~"(aligned-misaligned)")) +  # set the names for main, x and y axises
  theme(legend.position = "none") +
  NULL
# ggsave(filename = "ccf_d_fi.pdf", plot_ccf_d_fi, width = 4, height = 6)
plot_ccf_d_fi

2.4.2 Correct response times

2.4.2.1 Maximal model

# thisfile <- here("jubail", "lmm_ccf_rt_max.rds")
# 
# if (file.exists(thisfile)) {
#   lmm_ccf_rt_max <- readRDS(thisfile)
#   
# } else {
#   
#   lmm_ccf_rt_max <-
#     lmer(log(RT) ~ Congruency * Alignment + SD + Cue +
#            (Congruency * Alignment | Subject) +
#            (Congruency * Alignment | StimGroup),
#          filter(df_ccf, Correct),
#          control = lmerControl(optCtrl = list(maxfun = 1e7)))
#   
#   saveRDS(lmm_ccf_rt_max, thisfile)
# }
# 
# summary(lmm_ccf_rt_max)

2.4.2.2 Zero-correlation-parameter model

thisfile <- here("jubail", "lmm_ccf_rt_zcp.rds")

if (file.exists(thisfile)) {
  lmm_ccf_rt_zcp <- readRDS(thisfile)
  
} else {
  
  lmm_ccf_rt_zcp <- 
    lmer(log(RT) ~ Congruency * Alignment + SD + Cue + 
           (Con_C + Ali_C + Con_Ali || Subject) +
           (Con_C + Ali_C + Con_Ali || StimGroup),
         filter(df_ccf, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_ccf_rt_zcp, thisfile)
}

summary(lmm_ccf_rt_zcp)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ Congruency * Alignment + SD + Cue + (Con_C + Ali_C +      Con_Ali || Subject) + (Con_C + Ali_C + Con_Ali || StimGroup)
##    Data: filter(df_ccf, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 68232.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -5.0846 -0.6759 -0.1338  0.5492  9.4302 
## 
## Random effects:
##  Groups      Name        Variance  Std.Dev. 
##  Subject     (Intercept) 3.916e-02 1.979e-01
##  Subject.1   Con_C       9.291e-04 3.048e-02
##  Subject.2   Ali_C       1.189e-03 3.448e-02
##  Subject.3   Con_Ali     1.331e-03 3.648e-02
##  StimGroup   (Intercept) 4.113e-05 6.413e-03
##  StimGroup.1 Con_C       4.024e-05 6.343e-03
##  StimGroup.2 Ali_C       9.466e-05 9.729e-03
##  StimGroup.3 Con_Ali     1.438e-09 3.792e-05
##  Residual                1.063e-01 3.261e-01
## Number of obs: 109935, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                              Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)                 6.811e+00  9.549e-03  4.247e+02 713.236  < 2e-16 ***
## Congruency2-1               3.961e-02  3.160e-03  1.433e+01  12.535 4.07e-09 ***
## Alignment2-1                8.895e-03  4.000e-03  1.285e+01   2.224   0.0447 *  
## SD2-1                       7.212e-02  1.984e-03  1.086e+05  36.348  < 2e-16 ***
## Cue2-1                      4.191e-02  1.976e-03  1.084e+05  21.211  < 2e-16 ***
## Congruency2-1:Alignment2-1 -4.053e-02  4.308e-03  4.397e+02  -9.409  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Alg2-1 SD2-1  Cue2-1
## Congrncy2-1  0.004                            
## Alignmnt2-1  0.000 -0.009                     
## SD2-1        0.006 -0.005  0.011              
## Cue2-1       0.003  0.002  0.001  0.025       
## Cng2-1:A2-1 -0.003 -0.006  0.031  0.013  0.009
## optimizer (nloptwrap) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.00906044 (tol = 0.002, component 1)
summary(rePCA(lmm_ccf_rt_zcp))
## $Subject
## Importance of components:
##                          [,1]    [,2]   [,3]    [,4]
## Standard deviation     0.6069 0.11188 0.1057 0.09348
## Proportion of Variance 0.9191 0.03123 0.0279 0.02180
## Cumulative Proportion  0.9191 0.95029 0.9782 1.00000
## 
## $StimGroup
## Importance of components:
##                           [,1]    [,2]    [,3]      [,4]
## Standard deviation     0.02984 0.01967 0.01945 0.0001163
## Proportion of Variance 0.53775 0.23365 0.22859 0.0000100
## Cumulative Proportion  0.53775 0.77140 0.99999 1.0000000

By-item Con_Ali was removed.

2.4.2.3 Reduced model

thisfile <- here("jubail", "lmm_ccf_rt_rdc.rds")

if (file.exists(thisfile)) {
  lmm_ccf_rt_rdc <- readRDS(thisfile)
  
} else {
  
  lmm_ccf_rt_rdc <- 
    lmer(log(RT) ~ Congruency * Alignment + SD + Cue + 
           (Con_C + Ali_C + Con_Ali || Subject) +
           (Con_C + Ali_C || StimGroup), # + Con_Ali
         filter(df_ccf, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_ccf_rt_rdc, thisfile)
}

summary(lmm_ccf_rt_rdc)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ Congruency * Alignment + SD + Cue + (Con_C + Ali_C +      Con_Ali || Subject) + (Con_C + Ali_C || StimGroup)
##    Data: filter(df_ccf, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 68232.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -5.0846 -0.6759 -0.1338  0.5491  9.4301 
## 
## Random effects:
##  Groups      Name        Variance  Std.Dev.
##  Subject     (Intercept) 3.913e-02 0.197804
##  Subject.1   Con_C       9.304e-04 0.030503
##  Subject.2   Ali_C       1.189e-03 0.034484
##  Subject.3   Con_Ali     1.331e-03 0.036487
##  StimGroup   (Intercept) 4.110e-05 0.006411
##  StimGroup.1 Con_C       3.996e-05 0.006321
##  StimGroup.2 Ali_C       9.469e-05 0.009731
##  Residual                1.063e-01 0.326091
## Number of obs: 109935, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                              Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)                 6.811e+00  9.544e-03  4.252e+02 713.563  < 2e-16 ***
## Congruency2-1               3.961e-02  3.156e-03  1.441e+01  12.550 3.75e-09 ***
## Alignment2-1                8.895e-03  4.001e-03  1.284e+01   2.223   0.0448 *  
## SD2-1                       7.212e-02  1.984e-03  1.086e+05  36.348  < 2e-16 ***
## Cue2-1                      4.191e-02  1.976e-03  1.084e+05  21.211  < 2e-16 ***
## Congruency2-1:Alignment2-1 -4.053e-02  4.308e-03  4.404e+02  -9.408  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Alg2-1 SD2-1  Cue2-1
## Congrncy2-1  0.004                            
## Alignmnt2-1  0.000 -0.009                     
## SD2-1        0.006 -0.005  0.011              
## Cue2-1       0.003  0.002  0.001  0.025       
## Cng2-1:A2-1 -0.003 -0.006  0.031  0.013  0.009
## optimizer (nloptwrap) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.0228933 (tol = 0.002, component 1)
summary(rePCA(lmm_ccf_rt_rdc))
## $Subject
## Importance of components:
##                          [,1]    [,2]    [,3]    [,4]
## Standard deviation     0.6066 0.11189 0.10575 0.09354
## Proportion of Variance 0.9190 0.03127 0.02793 0.02185
## Cumulative Proportion  0.9190 0.95022 0.97815 1.00000
## 
## $StimGroup
## Importance of components:
##                           [,1]    [,2]    [,3]
## Standard deviation     0.02984 0.01966 0.01939
## Proportion of Variance 0.53880 0.23384 0.22736
## Cumulative Proportion  0.53880 0.77264 1.00000

2.4.2.4 Extended model

thisfile <- here("jubail", "lmm_ccf_rt_etd.rds")

if (file.exists(thisfile)) {
  lmm_ccf_rt_etd <- readRDS(thisfile)
  
} else {
  
  lmm_ccf_rt_etd <- 
    lmer(log(RT) ~ Congruency * Alignment + SD + Cue + 
           (Con_C + Ali_C + Con_Ali | Subject) +
           (Con_C + Ali_C | StimGroup), # + Con_Ali
         filter(df_ccf, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_ccf_rt_etd, thisfile)
}

summary(lmm_ccf_rt_etd)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ Congruency * Alignment + SD + Cue + (Con_C + Ali_C +      Con_Ali | Subject) + (Con_C + Ali_C | StimGroup)
##    Data: filter(df_ccf, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 68201.1
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -5.1052 -0.6756 -0.1345  0.5495  9.4508 
## 
## Random effects:
##  Groups    Name        Variance  Std.Dev. Corr             
##  Subject   (Intercept) 3.921e-02 0.198019                  
##            Con_C       9.519e-04 0.030853  0.14            
##            Ali_C       1.206e-03 0.034732  0.01 -0.08      
##            Con_Ali     1.420e-03 0.037685 -0.02 -0.97  0.29
##  StimGroup (Intercept) 4.103e-05 0.006405                  
##            Con_C       3.937e-05 0.006275 -0.39            
##            Ali_C       9.488e-05 0.009740  0.32  0.06      
##  Residual              1.063e-01 0.326083                  
## Number of obs: 109935, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                              Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)                 6.811e+00  9.554e-03  4.248e+02 712.859  < 2e-16 ***
## Congruency2-1               3.957e-02  3.154e-03  1.455e+01  12.544 3.38e-09 ***
## Alignment2-1                9.008e-03  4.008e-03  1.289e+01   2.248   0.0427 *  
## SD2-1                       7.219e-02  1.984e-03  1.086e+05  36.382  < 2e-16 ***
## Cue2-1                      4.193e-02  1.976e-03  1.079e+05  21.224  < 2e-16 ***
## Congruency2-1:Alignment2-1 -4.004e-02  4.330e-03  5.749e+02  -9.249  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Alg2-1 SD2-1  Cue2-1
## Congrncy2-1  0.014                            
## Alignmnt2-1  0.056  0.008                     
## SD2-1        0.006 -0.005  0.011              
## Cue2-1       0.003  0.002  0.001  0.025       
## Cng2-1:A2-1 -0.011 -0.189  0.080  0.013  0.008
## optimizer (nloptwrap) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.0158373 (tol = 0.002, component 1)
summary(rePCA(lmm_ccf_rt_etd))
## $Subject
## Importance of components:
##                          [,1]    [,2]    [,3]     [,4]
## Standard deviation     0.6074 0.15105 0.10315 0.003508
## Proportion of Variance 0.9168 0.05669 0.02644 0.000030
## Cumulative Proportion  0.9168 0.97353 0.99997 1.000000
## 
## $StimGroup
## Importance of components:
##                           [,1]    [,2]    [,3]
## Standard deviation     0.03091 0.02237 0.01389
## Proportion of Variance 0.57944 0.30356 0.11700
## Cumulative Proportion  0.57944 0.88300 1.00000

By-subject Con_C was removed.

thisfile <- here("jubail", "lmm_ccf_rt_etd2.rds")

if (file.exists(thisfile)) {
  lmm_ccf_rt_etd2 <- readRDS(thisfile)
  
} else {
  
  lmm_ccf_rt_etd2 <- 
    lmer(log(RT) ~ Congruency * Alignment + SD + Cue + 
           (Ali_C + Con_Ali | Subject) + # Con_C + 
           (Con_C + Ali_C | StimGroup), # + Con_Ali
         filter(df_ccf, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_ccf_rt_etd2, thisfile)
}

summary(lmm_ccf_rt_etd2)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ Congruency * Alignment + SD + Cue + (Ali_C + Con_Ali |      Subject) + (Con_C + Ali_C | StimGroup)
##    Data: filter(df_ccf, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 68276.2
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -5.0904 -0.6773 -0.1354  0.5492  9.4628 
## 
## Random effects:
##  Groups    Name        Variance  Std.Dev. Corr       
##  Subject   (Intercept) 3.915e-02 0.197859            
##            Ali_C       1.207e-03 0.034747  0.02      
##            Con_Ali     1.346e-03 0.036690 -0.01  0.28
##  StimGroup (Intercept) 4.066e-05 0.006376            
##            Con_C       9.661e-06 0.003108 -1.00      
##            Ali_C       9.456e-05 0.009724  0.29 -0.29
##  Residual              1.066e-01 0.326453            
## Number of obs: 109935, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                              Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)                 6.811e+00  9.545e-03  4.249e+02 713.549  < 2e-16 ***
## Congruency2-1               3.999e-02  2.208e-03  1.786e+01  18.110 6.11e-13 ***
## Alignment2-1                8.966e-03  4.005e-03  1.272e+01   2.239   0.0437 *  
## SD2-1                       7.212e-02  1.986e-03  1.085e+05  36.318  < 2e-16 ***
## Cue2-1                      4.199e-02  1.977e-03  1.078e+05  21.242  < 2e-16 ***
## Congruency2-1:Alignment2-1 -4.036e-02  4.315e-03  4.393e+02  -9.353  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Alg2-1 SD2-1  Cue2-1
## Congrncy2-1 -0.088                            
## Alignmnt2-1  0.054 -0.113                     
## SD2-1        0.006 -0.008  0.011              
## Cue2-1       0.003  0.002  0.001  0.025       
## Cng2-1:A2-1 -0.007 -0.009  0.077  0.013  0.008
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')
summary(rePCA(lmm_ccf_rt_etd2))
## $Subject
## Importance of components:
##                          [,1]    [,2]    [,3]
## Standard deviation     0.6061 0.12414 0.09244
## Proportion of Variance 0.9388 0.03938 0.02184
## Cumulative Proportion  0.9388 0.97816 1.00000
## 
## $StimGroup
## Importance of components:
##                           [,1]    [,2]    [,3]
## Standard deviation     0.03101 0.01995 2.7e-20
## Proportion of Variance 0.70717 0.29283 0.0e+00
## Cumulative Proportion  0.70717 1.00000 1.0e+00

By-item Con_C was removed.

thisfile <- here("jubail", "lmm_ccf_rt_etd3.rds")

if (file.exists(thisfile)) {
  lmm_ccf_rt_etd3 <- readRDS(thisfile)
  
} else {
  
  lmm_ccf_rt_etd3 <- 
    lmer(log(RT) ~ Congruency * Alignment + SD + Cue + 
           (Ali_C + Con_Ali | Subject) + # Con_C + 
           (Ali_C | StimGroup), # + Con_Ali + Con_C + 
         filter(df_ccf, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_ccf_rt_etd3, thisfile)
}

summary(lmm_ccf_rt_etd3)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ Congruency * Alignment + SD + Cue + (Ali_C + Con_Ali |      Subject) + (Ali_C | StimGroup)
##    Data: filter(df_ccf, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 68277.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -5.0803 -0.6771 -0.1350  0.5496  9.4577 
## 
## Random effects:
##  Groups    Name        Variance  Std.Dev. Corr       
##  Subject   (Intercept) 3.914e-02 0.19784             
##            Ali_C       1.208e-03 0.03476   0.02      
##            Con_Ali     1.346e-03 0.03669  -0.01  0.28
##  StimGroup (Intercept) 4.211e-05 0.00649             
##            Ali_C       9.410e-05 0.00970  0.31       
##  Residual              1.066e-01 0.32646             
## Number of obs: 109935, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                              Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)                 6.811e+00  9.551e-03  4.234e+02 713.049   <2e-16 ***
## Congruency2-1               4.006e-02  1.977e-03  1.088e+05  20.256   <2e-16 ***
## Alignment2-1                8.960e-03  3.999e-03  1.292e+01   2.240   0.0433 *  
## SD2-1                       7.212e-02  1.986e-03  1.089e+05  36.316   <2e-16 ***
## Cue2-1                      4.199e-02  1.977e-03  1.086e+05  21.236   <2e-16 ***
## Congruency2-1:Alignment2-1 -4.036e-02  4.315e-03  4.393e+02  -9.353   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Alg2-1 SD2-1  Cue2-1
## Congrncy2-1  0.007                            
## Alignmnt2-1  0.057 -0.014                     
## SD2-1        0.006 -0.009  0.011              
## Cue2-1       0.003  0.003  0.001  0.025       
## Cng2-1:A2-1 -0.006 -0.010  0.077  0.013  0.008
## optimizer (nloptwrap) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.00633072 (tol = 0.002, component 1)
summary(rePCA(lmm_ccf_rt_etd3))
## $Subject
## Importance of components:
##                          [,1]    [,2]    [,3]
## Standard deviation     0.6060 0.12417 0.09243
## Proportion of Variance 0.9387 0.03941 0.02184
## Cumulative Proportion  0.9387 0.97816 1.00000
## 
## $StimGroup
## Importance of components:
##                           [,1]   [,2]
## Standard deviation     0.03071 0.0183
## Proportion of Variance 0.73806 0.2619
## Cumulative Proportion  0.73806 1.0000

2.4.2.5 Optimal model

anova(lmm_ccf_rt_rdc, lmm_ccf_rt_etd3, refit=FALSE)

lmm_ccf_rt_rdc is used as the optimal model.

lmm_ccf_rt_opt <- lmm_ccf_rt_rdc

2.4.2.6 Effects of interest

RT in each condition:

emm_ccf_rt <- emmeans(lmm_ccf_rt_opt, ~ Congruency + Alignment)
summary(emm_ccf_rt, type = "response")

Plot of RT:

emmip(lmm_ccf_rt_opt, Congruency ~ Alignment, CIs = TRUE, type = "response") 

Congruency effect of RT in the aligned condition (differences of log RT):

emm_con_rt <- contrast(emm_ccf_rt, "pairwise", by = "Alignment")
summary(emm_con_rt[1], side="<", infer=TRUE)

Congruency effect of RT in the aligned condition (ratio of RT):

summary(emm_con_rt[1], side="<", infer=TRUE, type="response") # exp(-0.0599)

Composite effect of RT (differences of log RT):

emm_ccfe_rt <- contrast(emm_ccf_rt, interaction="pairwise")
summary(emm_ccfe_rt, side="<", infer=TRUE)

Composite effect of RT (ratio of RT):

summary(emm_ccfe_rt, side="<", infer=TRUE, type = "response") # exp(-0.0405)
# facilitation and interference of d'
contrast(emm_ccf_rt, "pairwise", by="Congruency")
## Congruency = con:
##  contrast  estimate      SE  df z.ratio p.value
##  ali - mis  -0.0292 0.00448 Inf  -6.502  <.0001
## 
## Congruency = inc:
##  contrast  estimate      SE  df z.ratio p.value
##  ali - mis   0.0114 0.00460 Inf   2.471  0.0135
## 
## Results are averaged over the levels of: SD, Cue 
## Degrees-of-freedom method: asymptotic 
## Results are given on the log (not the response) scale.

2.4.2.7 Plot

plot_ccf_rt <- emm_ccf_rt %>% 
  summary(type = "response") %>% 
  as_tibble() %>% 
  mutate(Congruency = fct_recode(Congruency, congruent="con", incongruent="inc"),
         Alignment = fct_recode(Alignment, aligned="ali", misaligned="mis")) %>%
  ggplot(aes(y = response, x = Alignment, color = Congruency, group = Congruency)) +
  geom_point(size = 2) + # position = position_dodge(width = 0.1), 
  geom_line(aes(linetype = Congruency), # position = position_dodge(width = 0.1),
            size = 0.8) +
  scale_linetype_manual(values=c("solid", "dashed"))+
  scale_color_manual(values=two_colors) +
  geom_errorbar(aes(ymin = asymp.LCL, ymax = asymp.UCL), size=1.5, width=0, 
                alpha = .6, # position = position_dodge(width = 0.1),
                show.legend = F) + 
  coord_cartesian(ylim = c(800, 1000)) +  # set the limit for y axis c(0, 1100)
  labs(y = "Correct response times (ms)", fill = "Congruency") +  # set the names for main, x and y axises
  # geom_text(label = c("", "***", "", ""), 
  #           color = "red", 
  #           size = 6, nudge_y = 50, nudge_x = 0.5) + # add starts to the significant columns
  NULL
# ggsave(filename = "ccf_rt.pdf", plot_ccf_rt, width = 8, height = 6)
plot_ccf_rt

plot_ccf_rt_fi <- contrast(emm_ccf_rt, "pairwise", by="Congruency") %>% 
  summary(infer=TRUE) %>% 
  as_tibble() %>% 
  mutate(Congruency = fct_recode(Congruency, congruent="con", incongruent="inc")) %>% 
  ggplot(aes(y = estimate, x = Congruency, color = Congruency)) +
  geom_point(size = 2) +
  geom_errorbar(aes(ymin = asymp.LCL, ymax = asymp.UCL), size=1.5, width=0, 
                alpha = .6) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  scale_color_manual(values=two_colors) +
  # coord_cartesian(ylim = ylimit_cf_fi_d) +  # set the limit for y axis c(0, 1100)
  labs(x = "Congruency", y = "Correct response times (ms)") +  # set the names for main, x and y axises
  theme(legend.position = "none") +
  NULL
# ggsave(filename = "ccf_rt_fi.pdf", plot_ccf_rt_fi, width = 8, height = 6)
plot_ccf_rt_fi

2.5 Facilitation and interference in the complete composite face task

2.5.1 Sensitivity d’

2.5.1.1 Maximal model

# thisfile <- here("jubail", "lmm_fi_d_max.rds")
# 
# if (file.exists(thisfile)) {
#   lmm_fi_d_max <- readRDS(thisfile)
#   
# } else {
#   
#   lmm_fi_d_max <-
#     glmer(isSame ~ Congruency * SD + Cue +
#             (Congruency * SD | Subject) +
#             (Congruency * SD | StimGroup),
#           df_fi,
#           family = binomial(link = "probit"),
#           control = glmerControl(optCtrl = list(maxfun = 1e7)))
#   
#   saveRDS(lmm_fi_d_max, thisfile)
# }
# 
# summary(lmm_fi_d_max)

2.5.1.2 Zero-correlation-parameter model

thisfile <- here("jubail", "lmm_fi_d_zcp.rds")

if (file.exists(thisfile)) {
  lmm_fi_d_zcp <- readRDS(thisfile)
  
} else {
  
  lmm_fi_d_zcp <- 
    glmer(isSame ~ Congruency * SD + Cue +
            (ConInc + IncIso + SD_C + 
               ConInc_SD + IncIso_SD || Subject) +
            (ConInc + IncIso + SD_C + 
               ConInc_SD + IncIso_SD || StimGroup),
          df_fi,
          family = binomial(link = "probit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_fi_d_zcp, thisfile)
}

summary(lmm_fi_d_zcp)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( probit )
## Formula: isSame ~ Congruency * SD + Cue + (ConInc + IncIso + SD_C + ConInc_SD +      IncIso_SD || Subject) + (ConInc + IncIso + SD_C + ConInc_SD +      IncIso_SD || StimGroup)
##    Data: df_fi
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
## 106160.6 106342.6 -53061.3 106122.6   107116 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -7.2442 -0.5596  0.2067  0.5317  5.1292 
## 
## Random effects:
##  Groups      Name        Variance Std.Dev.
##  Subject     (Intercept) 0.067985 0.26074 
##  Subject.1   ConInc      0.002563 0.05062 
##  Subject.2   IncIso      0.051150 0.22616 
##  Subject.3   SD_C        0.305651 0.55286 
##  Subject.4   ConInc_SD   0.135539 0.36816 
##  Subject.5   IncIso_SD   0.032643 0.18067 
##  StimGroup   (Intercept) 0.008036 0.08964 
##  StimGroup.1 ConInc      0.001091 0.03303 
##  StimGroup.2 IncIso      0.001567 0.03958 
##  StimGroup.3 SD_C        0.030000 0.17320 
##  StimGroup.4 ConInc_SD   0.005055 0.07110 
##  StimGroup.5 IncIso_SD   0.005309 0.07286 
## Number of obs: 107135, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)          0.092904   0.031224   2.975  0.00293 ** 
## Congruency2-1       -0.193829   0.015582 -12.440  < 2e-16 ***
## Congruency3-2        0.014845   0.019589   0.758  0.44857    
## SD2-1               -1.670330   0.061357 -27.223  < 2e-16 ***
## Cue2-1               0.162148   0.008952  18.113  < 2e-16 ***
## Congruency2-1:SD2-1  1.115531   0.036324  30.710  < 2e-16 ***
## Congruency3-2:SD2-1 -0.663671   0.032648 -20.328  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Cng3-2 SD2-1  Cue2-1 C2-1:S
## Congrncy2-1 -0.017                                   
## Congrncy3-2  0.006 -0.169                            
## SD2-1       -0.002  0.010 -0.001                     
## Cue2-1      -0.001  0.002 -0.003 -0.008              
## Cn2-1:SD2-1  0.008 -0.051  0.001 -0.017  0.009       
## Cn3-2:SD2-1 -0.001  0.001 -0.006  0.009 -0.004 -0.174
summary(rePCA(lmm_fi_d_zcp))
## $Subject
## Importance of components:
##                          [,1]   [,2]   [,3]    [,4]    [,5]    [,6]
## Standard deviation     0.5529 0.3682 0.2607 0.22616 0.18067 0.05062
## Proportion of Variance 0.5132 0.2276 0.1142 0.08589 0.05481 0.00430
## Cumulative Proportion  0.5132 0.7408 0.8550 0.94088 0.99570 1.00000
## 
## $StimGroup
## Importance of components:
##                          [,1]    [,2]    [,3]    [,4]    [,5]    [,6]
## Standard deviation     0.1732 0.08964 0.07286 0.07110 0.03958 0.03303
## Proportion of Variance 0.5876 0.15739 0.10397 0.09901 0.03069 0.02137
## Cumulative Proportion  0.5876 0.74496 0.84893 0.94794 0.97863 1.00000

2.5.1.3 Extended model

thisfile <- here("jubail", "lmm_fi_d_etd.rds")

if (file.exists(thisfile)) {
  lmm_fi_d_etd <- readRDS(thisfile)
  
} else {
  
  lmm_fi_d_etd <- 
    glmer(isSame ~ Congruency * SD + Cue +
            (ConInc + IncIso + SD_C + 
               ConInc_SD + IncIso_SD | Subject) +
            (ConInc + IncIso + SD_C + 
               ConInc_SD + IncIso_SD | StimGroup),
          df_fi,
          family = binomial(link = "probit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_fi_d_etd, thisfile)
}

summary(lmm_fi_d_etd)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( probit )
## Formula: isSame ~ Congruency * SD + Cue + (ConInc + IncIso + SD_C + ConInc_SD +      IncIso_SD | Subject) + (ConInc + IncIso + SD_C + ConInc_SD +      IncIso_SD | StimGroup)
##    Data: df_fi
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
## 105893.1 106362.6 -52897.6 105795.1   107086 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -8.5956 -0.5622  0.1812  0.5391  5.6529 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev. Corr                         
##  Subject   (Intercept) 0.069239 0.26313                               
##            ConInc      0.011667 0.10801  -0.20                        
##            IncIso      0.057558 0.23991   0.09 -0.33                  
##            SD_C        0.326161 0.57111  -0.21  0.85 -0.14            
##            ConInc_SD   0.274164 0.52361  -0.03 -0.79 -0.05 -0.47      
##            IncIso_SD   0.113887 0.33747   0.03  0.89 -0.02  0.66 -0.97
##  StimGroup (Intercept) 0.007418 0.08613                               
##            ConInc      0.002799 0.05291   1.00                        
##            IncIso      0.002594 0.05093  -0.84 -0.84                  
##            SD_C        0.030910 0.17581   0.97  0.97 -0.90            
##            ConInc_SD   0.010425 0.10210  -0.77 -0.77  0.79 -0.68      
##            IncIso_SD   0.010573 0.10282   0.84  0.84 -0.92  0.81 -0.96
## Number of obs: 107135, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)          0.100306   0.030282   3.312 0.000925 ***
## Congruency2-1       -0.214827   0.021121 -10.171  < 2e-16 ***
## Congruency3-2        0.015484   0.022510   0.688 0.491539    
## SD2-1               -1.681609   0.062438 -26.932  < 2e-16 ***
## Cue2-1               0.162371   0.008963  18.115  < 2e-16 ***
## Congruency2-1:SD2-1  1.146368   0.047010  24.385  < 2e-16 ***
## Congruency3-2:SD2-1 -0.679294   0.042303 -16.058  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Cng3-2 SD2-1  Cue2-1 C2-1:S
## Congrncy2-1  0.674                                   
## Congrncy3-2 -0.517 -0.628                            
## SD2-1        0.737  0.783 -0.604                     
## Cue2-1      -0.001  0.001 -0.003 -0.007              
## Cn2-1:SD2-1 -0.469 -0.565  0.378 -0.537  0.007       
## Cn3-2:SD2-1  0.581  0.591 -0.519  0.666 -0.003 -0.806
## optimizer (Nelder_Mead) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')
summary(rePCA(lmm_fi_d_etd))
## $Subject
## Importance of components:
##                          [,1]   [,2]    [,3]    [,4]    [,5]     [,6]
## Standard deviation     0.7453 0.4232 0.25026 0.23523 0.01282 0.001228
## Proportion of Variance 0.6514 0.2100 0.07345 0.06489 0.00019 0.000000
## Cumulative Proportion  0.6514 0.8615 0.93491 0.99981 1.00000 1.000000
## 
## $StimGroup
## Importance of components:
##                          [,1]    [,2]    [,3]      [,4]      [,5]      [,6]
## Standard deviation     0.2396 0.08028 0.02925 0.0005061 0.0002439 6.077e-05
## Proportion of Variance 0.8872 0.09959 0.01322 0.0000000 0.0000000 0.000e+00
## Cumulative Proportion  0.8872 0.98678 1.00000 1.0000000 1.0000000 1.000e+00

By-subject ConInc, IncIso, by-item (Intercept), ConInc, and IncIso were removed.

thisfile <- here("jubail", "lmm_fi_d_etd2.rds")

if (file.exists(thisfile)) {
  lmm_fi_d_etd2 <- readRDS(thisfile)
  
} else {
  
  lmm_fi_d_etd2 <-
    glmer(isSame ~ Congruency * SD + Cue +
            (SD_C + # ConInc + IncIso +
               ConInc_SD + IncIso_SD | Subject) +
            (0 + SD_C + # ConInc + IncIso +
               ConInc_SD + IncIso_SD | StimGroup),
          df_fi,
          family = binomial(link = "probit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_fi_d_etd2, thisfile)
}

summary(lmm_fi_d_etd2)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( probit )
## Formula: isSame ~ Congruency * SD + Cue + (SD_C + ConInc_SD + IncIso_SD |      Subject) + (0 + SD_C + ConInc_SD + IncIso_SD | StimGroup)
##    Data: df_fi
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
## 106536.5 106756.8 -53245.2 106490.5   107112 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -7.6110 -0.5711  0.1962  0.5369  5.2672 
## 
## Random effects:
##  Groups    Name        Variance Std.Dev. Corr             
##  Subject   (Intercept) 0.067063 0.25897                   
##            SD_C        0.311743 0.55834  -0.19            
##            ConInc_SD   0.250981 0.50098  -0.09 -0.44      
##            IncIso_SD   0.110616 0.33259   0.06  0.64 -0.97
##  StimGroup SD_C        0.034726 0.18635                   
##            ConInc_SD   0.013542 0.11637  -0.78            
##            IncIso_SD   0.009653 0.09825   0.85 -0.99      
## Number of obs: 107135, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)          0.097922   0.013012   7.526 5.24e-14 ***
## Congruency2-1       -0.193740   0.011358 -17.057  < 2e-16 ***
## Congruency3-2        0.013925   0.010672   1.305    0.192    
## SD2-1               -1.659433   0.065185 -25.457  < 2e-16 ***
## Cue2-1               0.161831   0.008916  18.150  < 2e-16 ***
## Congruency2-1:SD2-1  1.129787   0.049256  22.937  < 2e-16 ***
## Congruency3-2:SD2-1 -0.658134   0.040815 -16.125  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Cng3-2 SD2-1  Cue2-1 C2-1:S
## Congrncy2-1 -0.061                                   
## Congrncy3-2  0.029 -0.429                            
## SD2-1       -0.078  0.013 -0.001                     
## Cue2-1      -0.003  0.004 -0.006 -0.007              
## Cn2-1:SD2-1 -0.023 -0.056  0.005 -0.622  0.006       
## Cn3-2:SD2-1  0.021  0.005 -0.013  0.687 -0.003 -0.844
## optimizer (Nelder_Mead) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.0281041 (tol = 0.002, component 1)
summary(rePCA(lmm_fi_d_etd2))
## $Subject
## Importance of components:
##                          [,1]   [,2]    [,3]    [,4]
## Standard deviation     0.7104 0.4190 0.24336 0.03079
## Proportion of Variance 0.6816 0.2371 0.07999 0.00128
## Cumulative Proportion  0.6816 0.9187 0.99872 1.00000
## 
## $StimGroup
## Importance of components:
##                          [,1]    [,2]      [,3]
## Standard deviation     0.2292 0.07333 0.0007496
## Proportion of Variance 0.9072 0.09284 0.0000100
## Cumulative Proportion  0.9072 0.99999 1.0000000

By-subject (Intercept) and by-item IncIso_SD were removed.

thisfile <- here("jubail", "lmm_fi_d_etd3.rds")

if (file.exists(thisfile)) {
  
  lmm_fi_d_etd3 <- readRDS(thisfile)
  
} else {
  lmm_fi_d_etd3 <-
    glmer(isSame ~ Congruency * SD + Cue +
            (0 + SD_C + # ConInc + IncIso +
               ConInc_SD + IncIso_SD | Subject) +
            (0 + SD_C + # ConInc + IncIso +
               ConInc_SD | StimGroup), # + IncIso_SD
          df_fi,
          family = binomial(link = "probit"),
          control = glmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_fi_d_etd3, thisfile)
}

summary(lmm_fi_d_etd3)
## Generalized linear mixed model fit by maximum likelihood (Laplace Approximation) ['glmerMod']
##  Family: binomial  ( probit )
## Formula: isSame ~ Congruency * SD + Cue + (0 + SD_C + ConInc_SD + IncIso_SD |      Subject) + (0 + SD_C + ConInc_SD | StimGroup)
##    Data: df_fi
## Control: glmerControl(optCtrl = list(maxfun = 1e+07))
## 
##      AIC      BIC   logLik deviance df.resid 
## 108822.1 108975.5 -54395.1 108790.1   107119 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -7.3751 -0.5810  0.2395  0.5559  4.1981 
## 
## Random effects:
##  Groups    Name      Variance Std.Dev. Corr       
##  Subject   SD_C      0.291470 0.5399              
##            ConInc_SD 0.249493 0.4995   -0.44      
##            IncIso_SD 0.098667 0.3141    0.65 -0.97
##  StimGroup SD_C      0.031953 0.1788              
##            ConInc_SD 0.004135 0.0643   -0.76      
## Number of obs: 107135, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                      Estimate Std. Error z value Pr(>|z|)    
## (Intercept)          0.088076   0.004452  19.784   <2e-16 ***
## Congruency2-1       -0.187221   0.010949 -17.099   <2e-16 ***
## Congruency3-2        0.013211   0.010383   1.272    0.203    
## SD2-1               -1.603505   0.062578 -25.624   <2e-16 ***
## Cue2-1               0.157460   0.008783  17.928   <2e-16 ***
## Congruency2-1:SD2-1  1.090084   0.038240  28.506   <2e-16 ***
## Congruency3-2:SD2-1 -0.631907   0.025660 -24.626   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Cng3-2 SD2-1  Cue2-1 C2-1:S
## Congrncy2-1 -0.143                                   
## Congrncy3-2  0.065 -0.435                            
## SD2-1       -0.011  0.012 -0.001                     
## Cue2-1      -0.009  0.004 -0.006 -0.007              
## Cn2-1:SD2-1  0.046 -0.060  0.002 -0.486  0.007       
## Cn3-2:SD2-1 -0.005  0.002 -0.012  0.160 -0.005 -0.544
summary(rePCA(lmm_fi_d_etd3))
## $Subject
## Importance of components:
##                          [,1]   [,2]    [,3]
## Standard deviation     0.6928 0.3994 0.01394
## Proportion of Variance 0.7503 0.2494 0.00030
## Cumulative Proportion  0.7503 0.9997 1.00000
## 
## $StimGroup
## Importance of components:
##                          [,1]    [,2]
## Standard deviation     0.1856 0.04054
## Proportion of Variance 0.9545 0.04553
## Cumulative Proportion  0.9545 1.00000

2.5.1.4 Optimal model

anova(lmm_fi_d_zcp, lmm_fi_d_etd3, refit=FALSE)

lmm_fi_d_zcp is used as the optimal model.

lmm_fi_d_opt <- lmm_fi_d_zcp

2.5.1.5 Effects of interest

Hit and false alarm in each condition:

(emm_fi_probit <- emmeans(lmm_fi_d_opt, ~ Congruency + SD, type = "response"))
##  Congruency SD    prob      SE  df asymp.LCL asymp.UCL
##  con        same 0.906 0.00798 Inf     0.889     0.920
##  inc        same 0.713 0.01546 Inf     0.682     0.742
##  iso        same 0.818 0.01245 Inf     0.793     0.842
##  con        diff 0.190 0.01267 Inf     0.166     0.215
##  inc        diff 0.303 0.01584 Inf     0.273     0.335
##  iso        diff 0.203 0.01329 Inf     0.178     0.230
## 
## Results are averaged over the levels of: Cue 
## Confidence level used: 0.95 
## Intervals are back-transformed from the probit scale

Sensitivity d’ in each condition:

emm_fi_d <- contrast(emm_fi_probit, "pairwise", simple="SD")
summary(emm_fi_d[1:3], infer=TRUE)

Plot of d’:

emmip(emm_fi_d, ~ Congruency, CIs = TRUE)

Facilitation in d’:

emm_fie_d <- contrast(emm_fi_probit, interaction="pairwise")
summary(emm_fie_d[2], side=">", infer=TRUE)

Interference in d’:

summary(emm_fie_d[3], side="<", infer=TRUE)

2.5.1.6 Plot

fi_colors <- c("con" = two_colors[1], "inc" = two_colors[2], "iso" = "black")
fi_linetype <- c("con" = "solid", "inc" = "dashed", "iso" = "solid")

plot_fi_d <- emm_fi_d %>%
  summary(infer=TRUE) %>% 
  as_tibble() %>% 
  mutate(Alignment = if_else(Congruency=="iso", "iso", "ali"),
         group_ai = if_else(Congruency=="inc", "NULL", "ai"),
         group_ii = if_else(Congruency=="con", "NULL", "ii")) %>%
  mutate(Alignment = fct_recode(Alignment, aligned="ali", isolated="iso")) %>%
  ggplot(aes(Alignment, estimate, color=Congruency, linetype=Congruency)) +
  geom_point(size = 2) + # position = position_dodge(width = 0.1),
  geom_line(aes(group=group_ai, color="con", linetype="con"), size = 0.8) +
  geom_line(aes(group=group_ii, color="inc", linetype="inc"), size = 0.8) +
  geom_errorbar(aes(ymin = asymp.LCL, ymax = asymp.UCL), size=1.5, width=0, 
                alpha = .6, # position = position_dodge(width = 0.1),
                linetype = "solid",
                show.legend = F) + 
  scale_color_manual(values=fi_colors, 
                     breaks = c("con", "inc"),
                     labels = c("Facilitation", "Interference", NULL)) +
  scale_linetype_manual(values = fi_linetype,
                     breaks = c("con", "inc"),
                     labels = c("Facilitation", "Interference", NULL)) +
  coord_cartesian(ylim = c(0,3.2)) +  # set the limit for y axis c(0, 1100)
  labs(y = expression("Sensitivity"~italic("d'")), color=NULL, linetype=NULL) +  # set the names for main, x and y axises
  # geom_text(label = c("***", "***", ""),
  #           color = "red",
  #           size = 6, nudge_y = c(0.5, -0.5, 0), nudge_x = 0.5) + # add starts to the significant columns
  NULL
# ggsave(filename = "fi_d.pdf", plot_fi_d, width = 8, height = 6)
plot_fi_d

2.5.2 Correct response times

2.5.2.1 Maximal model

# thisfile <- here("jubail", "lmm_fi_rt_max.rds")
# 
# if (file.exists(thisfile)) {
#   lmm_fi_rt_max <- readRDS(thisfile)
# 
# } else {
#   
#   lmm_fi_rt_max <-
#     lmer(log(RT) ~ Congruency + SD + Cue +
#            (Congruency | Subject) +
#            (Congruency | StimGroup),
#          filter(df_fi, Correct),
#          control = lmerControl(optCtrl = list(maxfun = 1e7)))
#   
#   saveRDS(lmm_fi_rt_max, thisfile)
# }
# 
# summary(lmm_fi_rt_max)

2.5.2.2 Zero-correlation-parameter model

thisfile <- here("jubail", "lmm_fi_rt_zcp.rds")

if (file.exists(thisfile)) {
  lmm_fi_rt_zcp <- readRDS(thisfile)
  
} else {
  
  lmm_fi_rt_zcp <- 
    lmer(log(RT) ~ Congruency + SD + Cue +
           (ConInc + IncIso || Subject) +
           (ConInc + IncIso || StimGroup),
         filter(df_fi, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_fi_rt_zcp, thisfile)
}

summary(lmm_fi_rt_zcp)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ Congruency + SD + Cue + (ConInc + IncIso || Subject) +      (ConInc + IncIso || StimGroup)
##    Data: filter(df_fi, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 45749
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -4.9619 -0.6659 -0.1246  0.5299  8.2400 
## 
## Random effects:
##  Groups      Name        Variance  Std.Dev.
##  Subject     (Intercept) 3.639e-02 0.190749
##  Subject.1   ConInc      1.335e-03 0.036534
##  Subject.2   IncIso      2.088e-03 0.045696
##  StimGroup   (Intercept) 4.143e-05 0.006437
##  StimGroup.1 ConInc      2.241e-06 0.001497
##  StimGroup.2 IncIso      4.372e-05 0.006612
##  Residual                9.889e-02 0.314473
## Number of obs: 82627, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                 Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)    6.783e+00  9.238e-03  4.114e+02  734.29  < 2e-16 ***
## Congruency2-1  6.028e-02  3.244e-03  2.324e+01   18.58  1.9e-15 ***
## Congruency3-2 -9.796e-02  4.069e-03  2.295e+01  -24.07  < 2e-16 ***
## SD2-1          6.319e-02  2.205e-03  8.154e+04   28.66  < 2e-16 ***
## Cue2-1         4.360e-02  2.197e-03  8.136e+04   19.84  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Cng3-2 SD2-1 
## Congrncy2-1  0.008                     
## Congrncy3-2 -0.004 -0.305              
## SD2-1        0.003 -0.016  0.001       
## Cue2-1       0.004 -0.004  0.002  0.032
summary(rePCA(lmm_fi_rt_zcp))
## $Subject
## Importance of components:
##                          [,1]    [,2]    [,3]
## Standard deviation     0.6066 0.14531 0.11617
## Proportion of Variance 0.9140 0.05245 0.03353
## Cumulative Proportion  0.9140 0.96647 1.00000
## 
## $StimGroup
## Importance of components:
##                           [,1]    [,2]     [,3]
## Standard deviation     0.02103 0.02047 0.004761
## Proportion of Variance 0.50028 0.47408 0.025650
## Cumulative Proportion  0.50028 0.97435 1.000000

No random effects need to be removed. Therefore, we apply the extended model (i.e., the maximal model in this case) directly.

2.5.2.3 Extended model

thisfile <- here("jubail", "lmm_fi_rt_etd.rds")

if (file.exists(thisfile)) {
  lmm_fi_rt_etd <- readRDS(thisfile)
  
} else {
  lmm_fi_rt_etd <- 
    lmer(log(RT) ~ Congruency + SD + Cue +
           (ConInc + IncIso | Subject) +
           (ConInc + IncIso | StimGroup),
         filter(df_fi, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_fi_rt_etd, thisfile)
}

summary(lmm_fi_rt_etd)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ Congruency + SD + Cue + (ConInc + IncIso | Subject) +      (ConInc + IncIso | StimGroup)
##    Data: filter(df_fi, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 45627.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -5.0082 -0.6641 -0.1250  0.5304  8.2651 
## 
## Random effects:
##  Groups    Name        Variance  Std.Dev. Corr       
##  Subject   (Intercept) 3.647e-02 0.190977            
##            ConInc      2.664e-03 0.051613  0.07      
##            IncIso      3.682e-03 0.060675 -0.40 -0.74
##  StimGroup (Intercept) 4.033e-05 0.006351            
##            ConInc      5.533e-05 0.007439 -0.61      
##            IncIso      1.079e-04 0.010389  0.22 -0.91
##  Residual              9.873e-02 0.314219            
## Number of obs: 82627, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                 Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)    6.783e+00  9.242e-03  4.137e+02  733.94  < 2e-16 ***
## Congruency2-1  5.946e-02  4.331e-03  2.321e+01   13.73 1.25e-12 ***
## Congruency3-2 -9.748e-02  5.146e-03  1.850e+01  -18.94 1.45e-13 ***
## SD2-1          6.313e-02  2.203e-03  8.149e+04   28.66  < 2e-16 ***
## Cue2-1         4.354e-02  2.196e-03  8.105e+04   19.83  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Cng3-2 SD2-1 
## Congrncy2-1 -0.026                     
## Congrncy3-2 -0.185 -0.724              
## SD2-1        0.003 -0.012  0.001       
## Cue2-1       0.004 -0.003  0.002  0.032
## optimizer (nloptwrap) convergence code: 0 (OK)
## boundary (singular) fit: see help('isSingular')
summary(rePCA(lmm_fi_rt_etd))
## $Subject
## Importance of components:
##                          [,1]   [,2]    [,3]
## Standard deviation     0.6133 0.2258 0.08104
## Proportion of Variance 0.8673 0.1176 0.01514
## Cumulative Proportion  0.8673 0.9849 1.00000
## 
## $StimGroup
## Importance of components:
##                           [,1]   [,2]      [,3]
## Standard deviation     0.04067 0.0202 1.672e-05
## Proportion of Variance 0.80206 0.1979 0.000e+00
## Cumulative Proportion  0.80206 1.0000 1.000e+00

By-item (Intercept) was removed.

thisfile <- here("jubail", "lmm_fi_rt_etd2.rds")

if (file.exists(thisfile)) {
  lmm_fi_rt_etd2 <- readRDS(thisfile)
  
} else {
  
  lmm_fi_rt_etd2 <- 
    lmer(log(RT) ~ Congruency + SD + Cue +
           (ConInc + IncIso | Subject) +
           (0 + ConInc + IncIso | StimGroup), 
         filter(df_fi, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_fi_rt_etd2, thisfile)
}

summary(lmm_fi_rt_etd2)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ Congruency + SD + Cue + (ConInc + IncIso | Subject) +      (0 + ConInc + IncIso | StimGroup)
##    Data: filter(df_fi, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 45647
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -4.9939 -0.6636 -0.1249  0.5312  8.2793 
## 
## Random effects:
##  Groups    Name        Variance  Std.Dev. Corr       
##  Subject   (Intercept) 3.647e-02 0.190964            
##            ConInc      2.661e-03 0.051585  0.07      
##            IncIso      3.675e-03 0.060620 -0.40 -0.73
##  StimGroup ConInc      5.501e-05 0.007417            
##            IncIso      1.116e-04 0.010565 -1.00      
##  Residual              9.877e-02 0.314282            
## Number of obs: 82627, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                 Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)    6.783e+00  9.021e-03  4.538e+02  751.94  < 2e-16 ***
## Congruency2-1  5.944e-02  4.327e-03  2.087e+01   13.73 6.38e-12 ***
## Congruency3-2 -9.745e-02  5.181e-03  1.832e+01  -18.81 1.98e-13 ***
## SD2-1          6.292e-02  2.202e-03  8.168e+04   28.57  < 2e-16 ***
## Cue2-1         4.342e-02  2.195e-03  8.154e+04   19.78  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Cng3-2 SD2-1 
## Congrncy2-1  0.046                     
## Congrncy3-2 -0.220 -0.755              
## SD2-1        0.003 -0.012  0.000       
## Cue2-1       0.004 -0.003  0.002  0.031
## optimizer (nloptwrap) convergence code: 0 (OK)
## Model failed to converge with max|grad| = 0.00430193 (tol = 0.002, component 1)
summary(rePCA(lmm_fi_rt_etd2))
## $Subject
## Importance of components:
##                          [,1]   [,2]    [,3]
## Standard deviation     0.6131 0.2255 0.08111
## Proportion of Variance 0.8675 0.1174 0.01518
## Cumulative Proportion  0.8675 0.9848 1.00000
## 
## $StimGroup
## Importance of components:
##                           [,1]      [,2]
## Standard deviation     0.04107 9.808e-05
## Proportion of Variance 0.99999 1.000e-05
## Cumulative Proportion  0.99999 1.000e+00

By-subject ConInc was removed.

thisfile <- here("jubail", "lmm_fi_rt_etd3.rds")

if (file.exists(thisfile)) {
  lmm_fi_rt_etd3 <- readRDS(thisfile)
  
} else {
  
  lmm_fi_rt_etd3 <- 
    lmer(log(RT) ~ Congruency + SD + Cue +
           (ConInc + IncIso | Subject) +
           (0 + IncIso | StimGroup), # ConInc + 
         filter(df_fi, Correct),
         control = lmerControl(optCtrl = list(maxfun = 1e7)))
  
  saveRDS(lmm_fi_rt_etd3, thisfile)
}

summary(lmm_fi_rt_etd3)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
## Formula: log(RT) ~ Congruency + SD + Cue + (ConInc + IncIso | Subject) +      (0 + IncIso | StimGroup)
##    Data: filter(df_fi, Correct)
## Control: lmerControl(optCtrl = list(maxfun = 1e+07))
## 
## REML criterion at convergence: 45650.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -5.0001 -0.6636 -0.1243  0.5305  8.2827 
## 
## Random effects:
##  Groups    Name        Variance  Std.Dev. Corr       
##  Subject   (Intercept) 3.647e-02 0.190977            
##            ConInc      2.663e-03 0.051603  0.07      
##            IncIso      3.681e-03 0.060667 -0.40 -0.73
##  StimGroup IncIso      4.146e-05 0.006439            
##  Residual              9.878e-02 0.314294            
## Number of obs: 82627, groups:  Subject, 455; StimGroup, 10
## 
## Fixed effects:
##                 Estimate Std. Error         df t value Pr(>|t|)    
## (Intercept)    6.783e+00  9.022e-03  4.538e+02  751.89   <2e-16 ***
## Congruency2-1  5.958e-02  3.637e-03  4.399e+02   16.38   <2e-16 ***
## Congruency3-2 -9.757e-02  4.454e-03  3.701e+01  -21.91   <2e-16 ***
## SD2-1          6.297e-02  2.202e-03  8.167e+04   28.59   <2e-16 ***
## Cue2-1         4.349e-02  2.195e-03  8.155e+04   19.81   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Cng2-1 Cng3-2 SD2-1 
## Congrncy2-1  0.055                     
## Congrncy3-2 -0.256 -0.562              
## SD2-1        0.003 -0.014  0.001       
## Cue2-1       0.004 -0.004  0.002  0.031

2.5.2.4 Optimal model

anova(lmm_fi_rt_zcp, lmm_fi_rt_etd3, refit=FALSE)

lmm_fi_rt_etd3 is used as the optimal model.

lmm_fi_rt_opt <- lmm_fi_rt_etd3

2.5.2.5 Effects of interest

Correct RT for each condition:

emm_fi_rt <- emmeans(lmm_fi_rt_opt, ~ Congruency)
summary(emm_fi_rt, type = "response")

Plot of each condition:

emmip(emm_fi_rt,  ~ Congruency, CIs = TRUE, type = "response") 

Facilitation in RT (differences of log RT):

emm_fie_rt <- contrast(emm_fi_rt, "pairwise")
summary(emm_fie_rt[2], side="<", infer=TRUE)

Facilitation in RT (ratio of RT):

summary(emm_fie_rt[2], side="<", infer=TRUE, type="response") # exp(0.038)

Interference in RT (differences of log RT):

summary(emm_fie_rt[3], side=">", infer=TRUE) 

Interference in RT (ratio of RT):

summary(emm_fie_rt[3], side=">", infer=TRUE, type="response") # exp(0.0976)

2.5.2.6 Plot

fi_colors <- c("con" = two_colors[1], "inc" = two_colors[2], "iso" = "black")
fi_linetype <- c("con" = "solid", "inc" = "dashed", "iso" = "solid")

plot_fi_rt <- emm_fi_rt %>% 
  summary(type = "response") %>% 
  as_tibble() %>% 
  mutate(Alignment = if_else(Congruency=="iso", "iso", "ali"),
         group_ai = if_else(Congruency=="inc", "NULL", "ai"),
         group_ii = if_else(Congruency=="con", "NULL", "ii")) %>%
  mutate(Alignment = fct_recode(Alignment, aligned="ali", isolated="iso")) %>%
  ggplot(aes(y = response, x = Alignment, color=Congruency, linetype=Congruency)) +
  geom_point(size = 2) + # position = position_dodge(width = 0.1), color=c(two_colors, "black")
  geom_line(aes(group=group_ai, color="con", linetype="con"), size = 0.8) +
  geom_line(aes(group=group_ii, color="inc", linetype="inc"), size = 0.8) +
  geom_errorbar(aes(ymin = asymp.LCL, ymax = asymp.UCL), size=1.5, width=0, 
                alpha = .6, # position = position_dodge(width = 0.1),
                linetype = "solid",
                show.legend = F) + 
  scale_color_manual(values=fi_colors, 
                     breaks = c("con", "inc"),
                     labels = c("Facilitation", "Interference", NULL)) +
  scale_linetype_manual(values = fi_linetype,
                     breaks = c("con", "inc"),
                     labels = c("Facilitation", "Interference", NULL)) +
  coord_cartesian(ylim = c(800, 1000)) +  # set the limit for y axis c(0, 1100)
  labs(y = "Correct response times (ms)", color=NULL, linetype=NULL) +  # set the names for main, x and y axises
  # geom_text(label = c("***", "***", ""),
  #           color = "red",
  #           size = 6, nudge_y = c(-50, 20, 0), nudge_x = 0.5) + # add starts to the significant columns
  NULL
# ggsave(filename = "fi_rt.pdf", plot_fi_rt, width = 8, height = 6)
plot_fi_rt

2.6 Publication plot

2.6.1 Figure 2: holsitic processing effects

acc_y_limits <- c(.75, 1)
d_y_limits <- c(0, 2.7)
d_y_breaks <- seq(0, 2.5, .5)
rt_y_limits <- c(550, 1450)
rt_y_breaks <- seq(600, 1400, 200)

plot_pw_acc_pub <- plot_pw_acc + 
  coord_cartesian(ylim = acc_y_limits) 

plot_pw_rt_pub <- plot_pw_rt +
  scale_y_continuous(breaks=rt_y_breaks) +
  coord_cartesian(ylim = rt_y_limits) 
  
plot_scf_acc_pub <- plot_scf_acc +
  coord_cartesian(ylim = acc_y_limits) 

plot_scf_rt_pub <- plot_scf_rt +
  scale_y_continuous(breaks=rt_y_breaks) +
  coord_cartesian(ylim = rt_y_limits) 

plot_ccf_d_pub <- plot_ccf_d +
  coord_cartesian(ylim = d_y_limits) +
  scale_y_continuous(breaks=d_y_breaks) +
  theme(legend.position = c(.7, .2),
        legend.box = "horizontal") 

plot_ccf_rt_pub <- plot_ccf_rt +
  scale_y_continuous(breaks=rt_y_breaks) +
  coord_cartesian(ylim = rt_y_limits) +
  theme(legend.position = "none")

plot_fi_d_pub <- plot_fi_d +
  coord_cartesian(ylim = d_y_limits) +
  scale_y_continuous(breaks=d_y_breaks) +
  theme(legend.position = c(.7, .2),
        legend.box = "horizontal") 

plot_fi_rt_pub <- plot_fi_rt +
  scale_y_continuous(breaks=rt_y_breaks) +
  coord_cartesian(ylim = rt_y_limits) +
  theme(legend.position = "none")
plot_hp_iso <- 
  ggarrange(ggarrange(plot_pw_acc_pub, plot_pw_rt_pub) %>% 
              annotate_figure(top = text_grob("Part-whole task", 
                                              face="bold", size = 14)),
            ggplot() + theme_void(), 
            ggarrange(plot_scf_acc_pub, plot_scf_rt_pub) %>% 
              annotate_figure(top = text_grob("Standard composite face task", 
                                              face="bold", size = 14)),
            ggplot() + theme_void(), ggplot() + theme_void(), ggplot() + theme_void(), 
            ggarrange(plot_ccf_d_pub, plot_ccf_rt_pub) %>% 
              annotate_figure(top = text_grob("Complete composite face task", 
                                              face="bold", size = 14)),
            ggplot() + theme_void(), 
            ggarrange(plot_fi_d_pub, plot_fi_rt_pub) %>% 
              annotate_figure(top = text_grob("Facilitation and interference effects", 
                                              face="bold", size = 14)),
            labels = c("a", "", "b", "", "", "", "c", "", "d"),
            widths = c(1, 0.05, 1),
            heights = c(1, 0.05, 1),
            nrow=3, ncol=3) +
  bgcolor("white") +
  border("white")
# ggsave(filename = here("figures","plot_hp_iso.pdf"), plot_hp_iso, width = 14, height = 14*.65)
# plot_hp_iso
plot_hp_iso

plot_ccf_d_fi_pub <- plot_ccf_d_fi 
plot_ccf_rt_fi_pub <- plot_ccf_rt_fi 
plot_hp <- 
  ggarrange(ggarrange(plot_pw_acc_pub, plot_pw_rt_pub) %>% 
              annotate_figure(top = text_grob("Part-whole task", 
                                              face="bold", size = 14)),
            ggplot() + theme_void(), 
            ggarrange(plot_scf_acc_pub, plot_scf_rt_pub) %>% 
              annotate_figure(top = text_grob("Standard composite face task", 
                                              face="bold", size = 14)),
            ggplot() + theme_void(), ggplot() + theme_void(), ggplot() + theme_void(), 
            ggarrange(plot_ccf_d_pub, plot_ccf_rt_pub) %>% 
              annotate_figure(top = text_grob("Complete composite face task", 
                                              face="bold", size = 14)),
            ggplot() + theme_void(), 
            ggarrange(plot_ccf_d_fi_pub, plot_ccf_rt_fi_pub) %>% 
              annotate_figure(top = text_grob("Facilitation and interference effects", 
                                              face="bold", size = 14)),
            labels = c("a", "", "b", "", "", "", "c", "", "d"),
            widths = c(1, 0.05, 1),
            heights = c(1, 0.05, 1),
            nrow=3, ncol=3) +
  bgcolor("white") +
  border("white")
# ggsave(filename = here("figures","plot_hp.pdf"), plot_hp, width = 14, height = 14*.65)
# plot_hp
plot_hp

3 Relationships among effects

3.1 Calculate holistic processing effects

# custom function to use linear models to apply the regression method
cal_reg <- function(.data, formula, outvar){
  
  # .data: the input dataframe
  # formula: used in lm()
  # outvar: the column name for the regression effect (residuals) [new column]
  
  # perform linear regression
  lm_tmp <- lm(formula, .data)
  
  # save the residuals
  return(mutate(.data, {{outvar}} := residuals(lm_tmp))) # outvar is a variable
  # return(mutate(.data, "{outvar}" := residuals(lm_tmp))) # outvar is a string
}

3.1.1 Part-whole task

Condition means and part-whole effect in accuracy:

pwe_acc <- df_pw %>% 
  group_by(Subject, PW) %>% 
  summarize(acc = mean(Correct),
            .groups = "drop") %>% 
  pivot_wider(Subject, names_from = PW, values_from = acc) %>% 
  mutate(pwe_subt = whole - part) %>% # subtraction
  cal_reg(whole ~ part, pwe_regr) # regression

pwe_acc

Condition means and part-whole effect in RT:

pwe_rt <- df_pw %>% 
  filter(Correct) %>% 
  group_by(Subject, PW) %>% 
  summarize(rt = mean(RT),
            .groups = "drop") %>% 
  pivot_wider(Subject, names_from = PW, values_from = rt) %>% 
  mutate(pwe_subt = whole - part) %>% # subtraction
  cal_reg(whole ~ part, pwe_regr) # regression

pwe_rt

3.1.2 Standard composite face task

Condition means and composite effect in accuracy:

scfe_acc <- df_scf %>% 
  group_by(Subject, Alignment) %>% 
  summarize(acc = mean(Correct),
            .groups = "drop") %>% 
  pivot_wider(Subject, names_from = Alignment, values_from = acc) %>% 
  mutate(scfe_subt = ali - mis) %>% # subtraction
  cal_reg(ali ~ mis, scfe_regr) # regression

scfe_acc

Condition means and composite effect in RT:

scfe_rt <- df_scf %>% 
  filter(Correct) %>% 
  group_by(Subject, Alignment) %>% 
  summarize(rt = mean(RT),
            .groups = "drop") %>% 
  pivot_wider(Subject, names_from = Alignment, values_from = rt) %>% 
  mutate(scfe_subt = ali - mis) %>% # subtraction
  cal_reg(ali ~ mis, scfe_regr) # regression

scfe_rt

3.1.3 Complete composite face task

Condition means and composite effect in sensitivity d’:

ccfe_d <- df_ccf_iso %>% 
  group_by(Subject, Congruency, Alignment, SD) %>% 
  summarize(count = n(),
            same_adjust = (sum(isSame)+0.5)/(count+1), # Snodgrass & Corwin (1988)
            .groups = "drop") %>% 
  mutate(z = qnorm(same_adjust)) %>% 
  select(-c(same_adjust, count)) %>% 
  pivot_wider(names_from = "SD", values_from = z) %>% 
  mutate(d = same - diff) %>% 
  pivot_wider(Subject, names_from = c(Congruency, Alignment), values_from = d) %>% 
  rename(isolated = iso_iso) %>% 
  mutate(fac_subt = con_ali - con_mis, # subtraction
         int_subt = inc_ali - inc_mis, # subtraction
         fac_subt_iso = con_ali - isolated, # subtraction
         int_subt_iso = inc_ali - isolated, # subtraction
         cong_ali_subt = con_ali - inc_ali, # subtraction
         cong_mis_subt = con_mis - inc_mis, # subtraction
         ccfe_subt_subt = cong_ali_subt - cong_mis_subt) %>%  # subtraction
  cal_reg(con_ali ~ con_mis, fac_regr) %>% # regression
  cal_reg(inc_ali ~ inc_mis, int_regr) %>% # regression
  cal_reg(con_ali ~ isolated, fac_regr_iso) %>% # regression
  cal_reg(inc_ali ~ isolated, int_regr_iso) %>% # regression
  cal_reg(con_ali ~ inc_ali, cong_ali_regr) %>% # regression
  cal_reg(con_mis ~ inc_mis, cong_mis_regr) %>% # regression
  cal_reg(cong_ali_subt ~ cong_mis_subt, ccfe_subt_regr) %>% # regression
  cal_reg(cong_ali_regr ~ cong_mis_regr, ccfe_regr_regr) # regression

ccfe_d 

Condition means and composite effect in RT:

ccfe_rt <- df_ccf_iso %>% 
  filter(Correct) %>% 
  group_by(Subject, Congruency, Alignment) %>% 
  summarize(rt = mean(RT),
            .groups = "drop") %>% 
  pivot_wider(Subject, names_from = c(Congruency, Alignment), values_from = rt) %>% 
  rename(isolated = iso_iso) %>% 
  mutate(fac_subt = con_ali - con_mis, # subtraction
         int_subt = inc_ali - inc_mis, # subtraction
         fac_subt_iso = con_ali - isolated, # subtraction
         int_subt_iso = inc_ali - isolated, # subtraction
         cong_ali_subt = con_ali - inc_ali, # subtraction
         cong_mis_subt = con_mis - inc_mis, # subtraction
         ccfe_subt_subt = cong_ali_subt - cong_mis_subt) %>% # subtraction
  cal_reg(con_ali ~ con_mis, fac_regr) %>% # regression
  cal_reg(inc_ali ~ inc_mis, int_regr) %>% # regression
  cal_reg(con_ali ~ isolated, fac_regr_iso) %>% # regression
  cal_reg(inc_ali ~ isolated, int_regr_iso) %>% # regression
  cal_reg(con_ali ~ inc_ali, cong_ali_regr) %>% # regression
  cal_reg(con_mis ~ inc_mis, cong_mis_regr) %>% # regression
  cal_reg(cong_ali_subt ~ cong_mis_subt, ccfe_subt_regr) %>% # regression [first step is subtraction; second step is regression]
  cal_reg(cong_ali_regr ~ cong_mis_regr, ccfe_regr_regr) # regression

ccfe_rt

3.2 Reliability of measurement

psych::splitHalf(..., check.keys=FALSE) is used to calculate the Guttman’s \(\lambda_2\) (although it will report results of various methods).

# a general (custom) function to calculate the reliability for each condition
rel_each <- function(.data, DV, ...){
  
  out <- 
    .data %>%
    filter(...) %>% # to filter out conditions
    pivot_wider(Subject, values_from = all_of(DV), names_from = test_item) %>%
    select(-Subject) %>%
    data.matrix() %>% 
    splitHalf(check.keys=FALSE)
  
  return(out$lambda2) # only output lambda2
}
# function to calculate the reliability of subtraction scores
rel_sub <- function(re_inte, re_base, raw_inte, raw_base){
  
  # re_inte: the reliability of the condition of interest
  # re_base: the reliability of the baseline condition
  # raw_inte: the raw values for each participant in the condition of interest
  # raw_base: the raw values for each participant in the baseline condition
  
  tmp <- 2*sd(raw_inte)*sd(raw_base)*cor(raw_inte, raw_base)
  
  top <- var(raw_inte) * re_inte$Reliability + var(raw_base) * re_base$Reliability - tmp
  bot <- var(raw_inte) + var(raw_base) - tmp
  
  return(top/bot)
}
# function to calculate the reliability of regression scores
rel_reg <- function(re_inte, re_base, raw_inte, raw_base){
  
  # re_inte: the reliability of the condition of interest
  # re_base: the reliability of the baseline condition
  # raw_inte: the raw values for each participant in the condition of interest
  # raw_base: the raw values for each participant in the baseline condition
  
  tmp1 <- cor(raw_inte, raw_base)^2
  
  top <- re_inte$Reliability + re_base$Reliability * tmp1 - 2 * tmp1
  bot <- 1 - tmp1
  
  return(top/bot)
}
# function to calculate the upper boundary of the correlation
upper_boundary <- function(re1, re2){
  # reliability will be set to 0 if it is negative
  return(sqrt(max(re1$Reliability,0) * max(re2$Reliability,0)))
}

3.2.1 Part-whole task

Structure of part-whole task data for calculating reliability:

# PW
df_pw_rel <- df_pw %>%
  mutate(test_item = paste(str_remove(basename(StudyFace), ".png"), 
                           str_remove(basename(TestFace), ".png"), sep="_"))
str(df_pw_rel)
## tibble [42,880 × 12] (S3: tbl_df/tbl/data.frame)
##  $ Subject   : Factor w/ 480 levels "subj001","subj002",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ PW        : Factor w/ 2 levels "whole","part": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "whole" "part"
##   .. .. ..$ : chr "2-1"
##  $ Feature   : Factor w/ 2 levels "mouth","eyes": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "mouth" "eyes"
##   .. .. ..$ : chr "2-1"
##  $ Correct   : logi [1:42880] TRUE FALSE FALSE TRUE TRUE TRUE ...
##  $ RT        : num [1:42880] 2876 3001 1249 1439 1696 ...
##  $ StimGroup : chr [1:42880] "cau_m3" "cau_m5" "cau_m4" "cau_f5" ...
##  $ StudyFace : chr [1:42880] "stim/pw/main_stim/cau_m3.png" "stim/pw/main_stim/cau_m5.png" "stim/pw/main_stim/cau_m4.png" "stim/pw/main_stim/cau_f5.png" ...
##  $ TestFace  : chr [1:42880] "stim/pw/main_stim/cau_m3.mouth.whole.left.png" "stim/pw/main_stim/cau_m5.mouth.whole.left.png" "stim/pw/main_stim/cau_m4.mouth.whole.left.png" "stim/pw/main_stim/cau_f5.mouth.whole.right.png" ...
##  $ PW_C      : num [1:42880] -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 ...
##  $ Feature_C : num [1:42880] -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 -0.5 ...
##  $ PW_Feature: num [1:42880] 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 0.25 ...
##  $ test_item : chr [1:42880] "cau_m3_cau_m3.mouth.whole.left" "cau_m5_cau_m5.mouth.whole.left" "cau_m4_cau_m4.mouth.whole.left" "cau_f5_cau_f5.mouth.whole.right" ...

3.2.1.1 Reliability of each condition

Accuracy for whole

rel_pw_acc_w <- rel_each(df_pw_rel, "Correct", PW=="whole") %>%  
  tibble(DV = "acc", Condition = "whole", Approach = "GuttmanL2",  Reliability = .)
rel_pw_acc_w

Accuracy for part

rel_pw_acc_p <- rel_each(df_pw_rel, "Correct", PW=="part") %>% 
  tibble(DV = "acc", Condition = "part", Approach = "GuttmanL2", Reliability = .)
rel_pw_acc_p

RT for whole

rel_pw_rt_w <- rel_each(df_pw_rel, "RT", PW=="whole") %>% 
  tibble(DV = "rt", Condition = "whole", Approach = "GuttmanL2", Reliability = .)
rel_pw_rt_w

RT for part

rel_pw_rt_p <- rel_each(df_pw_rel, "RT", PW=="part") %>% 
  tibble(DV = "rt", Condition = "part", Approach = "GuttmanL2", Reliability = .)
rel_pw_rt_p

3.2.1.2 Reliability of the part-whole effect

Subtraction of accuracy

rel_pwe_acc_subt <- rel_sub(rel_pw_acc_w, rel_pw_acc_p,
                            pwe_acc$whole, pwe_acc$part) %>% 
  tibble(DV = "acc", Condition = "whole - part", 
         Approach = "Subtraction", Reliability = .)
rel_pwe_acc_subt

Regression of accuracy

rel_pwe_acc_regr <- rel_reg(rel_pw_acc_w, rel_pw_acc_p,
                            pwe_acc$whole, pwe_acc$part) %>% 
  tibble(DV = "acc", Condition = "whole ~ part", 
         Approach = "Regression", Reliability = .)
rel_pwe_acc_regr

Subtraction of RT

rel_pwe_rt_subt <- rel_sub(rel_pw_rt_w, rel_pw_rt_p,
                           pwe_rt$whole, pwe_rt$part) %>% 
  tibble(DV = "rt", Condition = "whole - part", 
         Approach = "Subtraction", Reliability = .)
rel_pwe_rt_subt

Regression of RT

rel_pwe_rt_regr <- rel_reg(rel_pw_rt_w, rel_pw_rt_p,
                           pwe_rt$whole, pwe_rt$part) %>% 
  tibble(DV = "rt", Condition = "whole ~ part", 
         Approach = "Regression", Reliability = .)
rel_pwe_rt_regr

3.2.1.3 Interim summary

rel_sum_pw <- bind_rows(rel_pw_acc_w, rel_pw_acc_p, rel_pw_rt_w, rel_pw_rt_p,
                        rel_pwe_acc_subt, rel_pwe_acc_regr, rel_pwe_rt_subt, rel_pwe_rt_regr) %>% 
  mutate(Task="PW", .before=1)
rel_sum_pw

3.2.2 Standard composite face task

Structure of standard composite face task data for calculating reliability:

# SCF
df_scf_rel <- df_scf %>% 
  group_by(Subject, Alignment) %>% 
  summarize(test_item = 1:n(),
            Correct, RT,
            .groups = "drop") 

str(df_scf_rel)
## tibble [35,648 × 5] (S3: tbl_df/tbl/data.frame)
##  $ Subject  : Factor w/ 480 levels "subj001","subj002",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ Alignment: Factor w/ 2 levels "ali","mis": 1 1 1 1 1 1 1 1 1 1 ...
##   ..- attr(*, "contrasts")= num [1:2, 1] -0.5 0.5
##   .. ..- attr(*, "dimnames")=List of 2
##   .. .. ..$ : chr [1:2] "ali" "mis"
##   .. .. ..$ : chr "2-1"
##  $ test_item: int [1:35648] 1 2 3 4 5 6 7 8 9 10 ...
##  $ Correct  : logi [1:35648] FALSE TRUE TRUE TRUE TRUE TRUE ...
##  $ RT       : num [1:35648] 749 528 817 798 1729 ...

3.2.2.1 Reliability of each condition

Accuracy for aligned

rel_scf_acc_a <- rel_each(df_scf_rel, "Correct", Alignment=="ali") %>% 
  tibble(DV = "acc", Condition = "aligned", Approach = "GuttmanL2",  Reliability = .)
rel_scf_acc_a

Accuracy for misaligned

rel_scf_acc_m <- rel_each(df_scf_rel, "Correct", Alignment=="mis") %>% 
  tibble(DV = "acc", Condition = "misaligned", Approach = "GuttmanL2",  Reliability = .)
rel_scf_acc_m

RT for aligned

rel_scf_rt_a <- rel_each(df_scf_rel, "RT", Alignment=="ali") %>% 
  tibble(DV = "rt", Condition = "aligned", Approach = "GuttmanL2",  Reliability = .)
rel_scf_rt_a

RT for misaligned

rel_scf_rt_m <- rel_each(df_scf_rel, "RT", Alignment=="mis") %>% 
  tibble(DV = "rt", Condition = "misaligned", Approach = "GuttmanL2",  Reliability = .)
rel_scf_rt_m

3.2.2.2 Reliability of the composite effect

Subtraction of accuracy

rel_scfe_acc_subt <- rel_sub(rel_scf_acc_a, rel_scf_acc_m, 
                             scfe_acc$ali, scfe_acc$mis) %>% 
  tibble(DV = "acc", Condition = "ali - mis", 
         Approach = "Subtraction", Reliability = .)
rel_scfe_acc_subt

Regression of accuracy

rel_scfe_acc_regr <- rel_reg(rel_scf_acc_a, rel_scf_acc_m, 
                             scfe_acc$ali, scfe_acc$mis) %>% 
  tibble(DV = "acc", Condition = "ali ~ mis", 
         Approach = "Regression", Reliability = .)
rel_scfe_acc_regr

Subtraction of RT

rel_scfe_rt_subt <- rel_sub(rel_scf_rt_a, rel_scf_rt_m, 
                            scfe_rt$ali, scfe_rt$mis) %>% 
  tibble(DV = "rt", Condition = "ali - mis", 
         Approach = "Subtraction", Reliability = .)
rel_scfe_rt_subt

Regression of RT

rel_scfe_rt_regr <- rel_reg(rel_scf_rt_a, rel_scf_rt_m, 
                            scfe_rt$ali, scfe_rt$mis) %>% 
  tibble(DV = "rt", Condition = "ali ~ mis", 
         Approach = "Regression", Reliability = .)
rel_scfe_rt_regr

3.2.2.3 Interim summary

rel_sum_scf <- bind_rows(rel_scf_acc_a, rel_scf_acc_m, rel_scf_rt_a, rel_scf_rt_m,
                         rel_scfe_acc_subt, rel_scfe_acc_regr, rel_scfe_rt_subt, rel_scfe_rt_regr) %>% 
  mutate(Task = "SCF", .before =1)
rel_sum_scf

3.2.3 Complete composite face task

Structure of complete composite face task data for calculating reliability:

# CCF
df_ccf_rel <- df_ccf_iso %>% # need to include the isolated condition
  group_by(Subject, Congruency, Alignment) %>% 
  summarize(test_item = 1:n(),
            Correct, RT,
            .groups = "drop") 

str(df_ccf_rel)
## tibble [178,640 × 6] (S3: tbl_df/tbl/data.frame)
##  $ Subject   : Factor w/ 480 levels "subj001","subj002",..: 1 1 1 1 1 1 1 1 1 1 ...
##  $ Congruency: Factor w/ 3 levels "con","inc","iso": 1 1 1 1 1 1 1 1 1 1 ...
##  $ Alignment : Factor w/ 3 levels "ali","mis","iso": 1 1 1 1 1 1 1 1 1 1 ...
##  $ test_item : int [1:178640] 1 2 3 4 5 6 7 8 9 10 ...
##  $ Correct   : logi [1:178640] TRUE TRUE TRUE TRUE TRUE TRUE ...
##  $ RT        : num [1:178640] 651 1399 986 819 1362 ...

3.2.3.1 Reliability of each condition (acc & RT)

Behavioral choices for congruent aligned

rel_ccf_beha_ca <- rel_each(df_ccf_rel, "Correct", 
                            Congruency=="con", Alignment=="ali") %>%
  tibble(DV = "acc", Condition = "con_ali", 
         Approach = "GuttmanL2",  Reliability = .)

rel_ccf_beha_ca

Behavioral choices for incongruent aligned

rel_ccf_beha_ia <- rel_each(df_ccf_rel, "Correct", 
                            Congruency=="inc", Alignment=="ali") %>% 
  tibble(DV = "acc", Condition = "inc_ali", 
         Approach = "GuttmanL2",  Reliability = .)
rel_ccf_beha_ia

Behavioral choices for congruent misaligned

rel_ccf_beha_cm <- rel_each(df_ccf_rel, "Correct", 
                            Congruency=="con", Alignment=="mis") %>%
  tibble(DV = "acc", Condition = "con_mis", 
         Approach = "GuttmanL2",  Reliability = .)
rel_ccf_beha_cm

Behavioral choices for incongruent misaligned

rel_ccf_beha_im <- rel_each(df_ccf_rel, "Correct", 
                            Congruency=="inc", Alignment=="mis") %>% 
  tibble(DV = "acc", Condition = "inc_mis", 
         Approach = "GuttmanL2",  Reliability = .)
rel_ccf_beha_im

Behavioral choices for isolated

rel_ccf_beha_iso <- rel_each(df_ccf_rel, "Correct", 
                             Congruency=="iso", Alignment=="iso") %>%
  tibble(DV = "acc", Condition = "isolated", 
         Approach = "GuttmanL2",  Reliability = .)
rel_ccf_beha_iso

RT for congruent aligned

rel_ccf_rt_ca <- rel_each(df_ccf_rel, "RT", 
                          Congruency=="con", Alignment=="ali") %>% 
  tibble(DV = "rt", Condition = "con_ali", 
         Approach = "GuttmanL2",  Reliability = .)
rel_ccf_rt_ca

RT for incongruent aligned

rel_ccf_rt_ia <- rel_each(df_ccf_rel, "RT", 
                          Congruency=="inc", Alignment=="ali") %>% 
  tibble(DV = "rt", Condition = "inc_ali", 
         Approach = "GuttmanL2",  Reliability = .)
rel_ccf_rt_ia

RT for congruent misaligned

rel_ccf_rt_cm <- rel_each(df_ccf_rel, "RT", 
                          Congruency=="con", Alignment=="mis") %>% 
  tibble(DV = "rt", Condition = "con_mis", 
         Approach = "GuttmanL2",  Reliability = .)
rel_ccf_rt_cm

RT for incongruent misaligned

rel_ccf_rt_im <- rel_each(df_ccf_rel, "RT", 
                          Congruency=="inc", Alignment=="mis") %>% 
  tibble(DV = "rt", Condition = "inc_mis", 
         Approach = "GuttmanL2",  Reliability = .)
rel_ccf_rt_im

RT for isolated

rel_ccf_rt_iso <- rel_each(df_ccf_rel, "RT", 
                           Congruency=="iso", Alignment=="iso") %>% 
  tibble(DV = "rt", Condition = "isolated", 
         Approach = "GuttmanL2",  Reliability = .)
rel_ccf_rt_iso

3.2.3.2 Reliability of each condition (d’)

rel_each_d <- function(.data, iter=100, startseed=12, cores=4){
  
  seeds <- startseed+1:iter
  ls_df <- pbapply::pblapply(seeds, rel_each_d_one, cl=cores,
                             df_each=.data)
  df_rel_d <- bind_rows(ls_df, .id = "iter") %>% 
    summarize(Reliability = mean(lambda2)) %>% 
    mutate(DV = "d",
           Approach = "GuttmanL2",
           .before=1) %>% 
    mutate(niter=iter)
  
  return(df_rel_d)
}

# calculate d' for one random split
rel_each_d_one <- function(df_each, seed){
  set.seed(seed)
  
  rel_value <- df_each %>% 
    # split data into two half
    group_by(Subject, Cue, SD) %>% 
    summarize(isSame,
              half = paste0("h", 2-(sample(n(), n())>n()/2)),
              .groups = "drop") %>% 
    # calculate d'
    group_by(Subject, SD, half) %>% 
    summarize(saysame = sum(isSame),
              count = n(),
              .groups = "drop") %>% 
    mutate(saysame = (saysame+.5)/(count+1), # correction
           saysame = qnorm(saysame)) %>% 
    select(-count) %>% 
    pivot_wider(names_from = SD, values_from = saysame) %>% 
    mutate(d = same - diff) %>% 
    select(-c(same, diff)) %>% 
    pivot_wider(names_from = half, values_from = d) %>% 
    select(-Subject) %>% 
    psych::splitHalf(check.keys=FALSE)
  
  return(data.frame(lambda2 = rel_value$lambda2))
}

# parameters for calculating d'
niter <- 10000
ncore <- 16

d’ for congruent aligned

thisfile <- file.path("simulation", "rel_ccf_d_ca.rds")

if (file.exists(thisfile)) {
  
  rel_ccf_d_ca <- readRDS(thisfile)
  
} else {
  rel_ccf_d_ca <- df_ccf_iso %>% 
    filter(Congruency=="con", Alignment=="ali") %>%
    rel_each_d(iter = niter, cores = ncore) %>% 
    mutate(Condition = "con_ali")
  
  saveRDS(rel_ccf_d_ca, thisfile)
}

rel_ccf_d_ca

d’ for incongruent aligned

thisfile <- file.path("simulation", "rel_ccf_d_ia.rds")

if (file.exists(thisfile)) {
  
  rel_ccf_d_ia <- readRDS(thisfile)
  
} else {
  rel_ccf_d_ia <- df_ccf_iso %>% 
    filter(Congruency=="inc", Alignment=="ali") %>%
    rel_each_d(iter = niter, cores = ncore) %>% 
    mutate(Condition = "inc_ali")
  
  saveRDS(rel_ccf_d_ia, thisfile)
}

rel_ccf_d_ia

d’ for congruent misaligned

thisfile <- file.path("simulation", "rel_ccf_d_cm.rds")

if (file.exists(thisfile)) {
  
  rel_ccf_d_cm <- readRDS(thisfile)
  
} else {
  rel_ccf_d_cm <- df_ccf_iso %>% 
    filter(Congruency=="con", Alignment=="mis") %>%
    rel_each_d(iter = niter, cores = ncore) %>% 
    mutate(Condition = "con_mis")
  
  saveRDS(rel_ccf_d_cm, thisfile)
}

rel_ccf_d_cm

d’ for incongruent misaligned

thisfile <- file.path("simulation", "rel_ccf_d_im.rds")

if (file.exists(thisfile)) {
  
  rel_ccf_d_im <- readRDS(thisfile)
  
} else {
  rel_ccf_d_im <- df_ccf_iso %>% 
    filter(Congruency=="inc", Alignment=="mis") %>%
    rel_each_d(iter = niter, cores = ncore) %>% 
    mutate(Condition="inc_mis")
  
  saveRDS(rel_ccf_d_im, thisfile)
}

rel_ccf_d_im

d’ for isolated

thisfile <- file.path("simulation", "rel_ccf_d_iso.rds")

if (file.exists(thisfile)) {
  
  rel_ccf_d_iso <- readRDS(thisfile)
  
} else {
  rel_ccf_d_iso <- df_ccf_iso %>% 
    filter(Congruency=="iso", Alignment=="iso") %>%
    rel_each_d(iter = niter, cores = ncore) %>% 
    mutate(Condition="isolated")
  
  saveRDS(rel_ccf_d_iso, thisfile)
}

rel_ccf_d_iso

3.2.3.3 Reliability of facilitation

Subtraction of accuracy (with isolated baselines)

rel_fac_d_iso_subt <- rel_sub(rel_ccf_d_ca, rel_ccf_d_iso, 
                              ccfe_d$con_ali, ccfe_d$isolated) %>% 
  tibble(DV = "d", Condition = "con_ali - isolated", 
         Approach = "Subtraction",  Reliability = .)

rel_fac_d_iso_subt

Regression of accuracy (with isolated baselines)

rel_fac_d_iso_regr <- rel_reg(rel_ccf_d_ca, rel_ccf_d_iso, 
                              ccfe_d$con_ali, ccfe_d$isolated) %>% 
  tibble(DV = "d", Condition = "con_ali ~ isolated", 
         Approach = "Regression",  Reliability = .)

rel_fac_d_iso_regr

Subtraction of RT (with isolated baselines)

rel_fac_rt_iso_subt <- rel_sub(rel_ccf_rt_ca, rel_ccf_rt_iso, 
                               ccfe_rt$con_ali, ccfe_rt$isolated) %>% 
  tibble(DV = "rt", Condition = "con_ali - isolated", 
         Approach = "Subtraction",  Reliability = .)

rel_fac_rt_iso_subt

Regression of RT (with isolated baselines)

rel_fac_rt_iso_regr <- rel_reg(rel_ccf_rt_ca, rel_ccf_rt_iso, 
                               ccfe_rt$con_ali, ccfe_rt$isolated) %>% 
  tibble(DV = "rt", Condition = "con_ali ~ isolated", 
         Approach = "Regression",  Reliability = .)

rel_fac_rt_iso_regr

Subtraction of accuracy (with misaligned baselines)

rel_fac_d_subt <- rel_sub(rel_ccf_d_ca, rel_ccf_d_cm, 
                          ccfe_d$con_ali, ccfe_d$con_mis) %>% 
  tibble(DV = "d", Condition = "con_ali - con_mis", 
         Approach = "Subtraction",  Reliability = .)

rel_fac_d_subt

Regression of accuracy (with misaligned baselines)

rel_fac_d_regr <- rel_reg(rel_ccf_d_ca, rel_ccf_d_cm, 
                          ccfe_d$con_ali, ccfe_d$con_mis) %>% 
  tibble(DV = "d", Condition = "con_ali ~ con_mis", 
         Approach = "Regression",  Reliability = .)

rel_fac_d_regr

Subtraction of RT (with misaligned baselines)

rel_fac_rt_subt <- rel_sub(rel_ccf_rt_ca, rel_ccf_rt_cm, 
                           ccfe_rt$con_ali, ccfe_rt$con_mis) %>% 
  tibble(DV = "rt", Condition = "con_ali - con_mis", 
         Approach = "Subtraction",  Reliability = .)

rel_fac_rt_subt

Regression of RT (with misaligned baselines)

rel_fac_rt_regr <- rel_reg(rel_ccf_rt_ca, rel_ccf_rt_cm, 
                           ccfe_rt$con_ali, ccfe_rt$con_mis) %>% 
  tibble(DV = "rt", Condition = "con_ali ~ con_mis", 
         Approach = "Regression",  Reliability = .)

rel_fac_rt_subt

3.2.3.4 Reliability of interference

Subtraction of d (with isolated baselines)

rel_int_d_iso_subt <- rel_sub(rel_ccf_d_ia, rel_ccf_d_iso, 
                              ccfe_d$inc_ali, ccfe_d$isolated) %>% 
  tibble(DV = "d", Condition = "inc_ali - isolated", 
         Approach = "Subtraction",  Reliability = .)

rel_int_d_iso_subt

Regression of d (with isolated baselines)

rel_int_d_iso_regr <- rel_reg(rel_ccf_d_ia, rel_ccf_d_iso,
                              ccfe_d$inc_ali, ccfe_d$isolated) %>% 
  tibble(DV = "d", Condition = "inc_ali ~ isolated", 
         Approach = "Regression",  Reliability = .)

rel_int_d_iso_regr

Subtraction of RT (with isolated baselines)

rel_int_rt_iso_subt <- rel_sub(rel_ccf_rt_ia, rel_ccf_rt_iso, 
                               ccfe_rt$inc_ali, ccfe_rt$isolated) %>% 
  tibble(DV = "rt", Condition = "inc_ali - isolated", 
         Approach = "Subtraction",  Reliability = .)

rel_int_rt_iso_subt

Regression of RT (with isolated baselines)

rel_int_rt_iso_regr <- rel_reg(rel_ccf_rt_ia, rel_ccf_rt_iso, 
                               ccfe_rt$inc_ali, ccfe_rt$isolated) %>% 
  tibble(DV = "rt", Condition = "inc_ali ~ isolated", 
         Approach = "Regression",  Reliability = .)

rel_int_rt_iso_regr

Subtraction of d (with misaligned baselines)

rel_int_d_subt <- rel_sub(rel_ccf_d_ia, rel_ccf_d_im, 
                          ccfe_d$inc_ali, ccfe_d$inc_mis) %>% 
  tibble(DV = "d", Condition = "inc_ali - inc_mis", 
         Approach = "Subtraction",  Reliability = .)

rel_int_d_subt

Regression of d (with misaligned baselines)

rel_int_d_regr <- rel_reg(rel_ccf_d_ia, rel_ccf_d_im,
                          ccfe_d$inc_ali, ccfe_d$inc_mis) %>% 
  tibble(DV = "d", Condition = "inc_ali ~ inc_mis", 
         Approach = "Regression",  Reliability = .)

rel_int_d_regr

Subtraction of RT (with misaligned baselines)

rel_int_rt_subt <- rel_sub(rel_ccf_rt_ia, rel_ccf_rt_im, 
                           ccfe_rt$inc_ali, ccfe_rt$inc_mis) %>% 
  tibble(DV = "rt", Condition = "inc_ali - inc_mis", 
         Approach = "Subtraction",  Reliability = .)

rel_int_rt_subt

Regression of RT (with misaligned baselines)

rel_int_rt_regr <- rel_reg(rel_ccf_rt_ia, rel_ccf_rt_im, 
                           ccfe_rt$inc_ali, ccfe_rt$inc_mis) %>% 
  tibble(DV = "rt", Condition = "inc_ali ~ inc_mis", 
         Approach = "Regression",  Reliability = .)

rel_int_rt_regr

3.2.3.5 Reliability of congruency effects

Congruency in d subtraction (aligned)

rel_ccf_d_cong_a_subt <- rel_sub(rel_ccf_d_ca, rel_ccf_d_ia, 
                                 ccfe_d$con_ali, ccfe_d$inc_ali) %>% 
  tibble(DV = "d", Condition = "con_ali - inc_ali", 
         Approach = "Subtraction",  Reliability = .)

rel_ccf_d_cong_a_subt

Congruency in d regression (aligned)

rel_ccf_d_cong_a_regr <- rel_reg(rel_ccf_d_ca, rel_ccf_d_ia, 
                                 ccfe_d$con_ali, ccfe_d$inc_ali) %>% 
  tibble(DV = "d", Condition = "con_ali ~ inc_ali", 
         Approach = "Regression",  Reliability = .)

rel_ccf_d_cong_a_regr

Congruency in d subtraction (misaligned)

rel_ccf_d_cong_m_subt <- rel_sub(rel_ccf_d_cm, rel_ccf_d_im, 
                                 ccfe_d$con_mis, ccfe_d$inc_mis) %>% 
  tibble(DV = "d", Condition = "con_mis - inc_mis", 
         Approach = "Subtraction",  Reliability = .)

rel_ccf_d_cong_m_subt

Congruency in d regression (misaligned)

rel_ccf_d_cong_m_regr <- rel_reg(rel_ccf_d_cm, rel_ccf_d_im, 
                                 ccfe_d$con_mis, ccfe_d$inc_mis) %>% 
  tibble(DV = "d", Condition = "con_mis ~ inc_mis", 
         Approach = "Regression",  Reliability = .)

rel_ccf_d_cong_m_regr

Congruency in RT subtraction (aligned)

rel_ccf_rt_cong_a_subt <- rel_sub(rel_ccf_rt_ca, rel_ccf_rt_ia, 
                                  ccfe_d$con_ali, ccfe_d$inc_ali) %>% 
  tibble(DV = "rt", Condition = "con_ali - inc_ali", 
         Approach = "Subtraction",  Reliability = .)

rel_ccf_rt_cong_a_subt

Congruency in RT regression (aligned)

rel_ccf_rt_cong_a_regr <- rel_reg(rel_ccf_rt_ca, rel_ccf_rt_ia, 
                                  ccfe_d$con_ali, ccfe_d$inc_ali) %>% 
  tibble(DV = "rt", Condition = "con_ali ~ inc_ali", 
         Approach = "Regression",  Reliability = .)

rel_ccf_rt_cong_a_regr

Congruency in RT subtraction (misaligned)

rel_ccf_rt_cong_m_subt <- rel_sub(rel_ccf_rt_cm, rel_ccf_rt_im, 
                                  ccfe_d$con_mis, ccfe_d$inc_mis) %>% 
  tibble(DV = "rt", Condition = "con_mis - inc_mis", 
         Approach = "Subtraction",  Reliability = .)

rel_ccf_rt_cong_m_subt

Congruency in RT subtraction (misaligned)

rel_ccf_rt_cong_m_regr <- rel_reg(rel_ccf_rt_cm, rel_ccf_rt_im, 
                                  ccfe_d$con_mis, ccfe_d$inc_mis) %>% 
  tibble(DV = "rt", Condition = "con_mis ~ inc_mis", 
         Approach = "Regression",  Reliability = .)

rel_ccf_rt_cong_m_regr

3.2.3.6 Reliability of the composite effect

Subtraction of Congruency (composite effect) in d [subtraction of subtraction]

rel_ccfe_d_subt_subt <- rel_sub(rel_ccf_d_cong_a_subt, rel_ccf_d_cong_m_subt, 
                                ccfe_d$cong_ali_subt, ccfe_d$cong_mis_subt) %>% 
  tibble(DV = "d", Condition = "(con_ali - inc_ali) - (con_mis - inc_mis)", 
         Approach = "Subtraction", Reliability = .)

rel_ccfe_d_subt_subt

Regression of Congruency (composite effect) in d [regression of subtraction]

rel_ccfe_d_subt_regr <- rel_reg(rel_ccf_d_cong_a_subt, rel_ccf_d_cong_m_subt, 
                                ccfe_d$cong_ali_subt, ccfe_d$cong_mis_subt) %>% 
  tibble(DV = "d", Condition = "(con_ali - inc_ali) ~ (con_mis - inc_mis)", 
         Approach = "Regression", Reliability = .)

rel_ccfe_d_subt_regr

Subtraction of Congruency (composite effect) in d [subtraction of regression]

rel_ccfe_d_regr_subt <- rel_sub(rel_ccf_d_cong_a_regr, rel_ccf_d_cong_m_regr, 
                                ccfe_d$cong_ali_regr, ccfe_d$cong_mis_regr) %>% 
  tibble(DV = "d", Condition = "(con_ali ~ inc_ali) - (con_mis ~ inc_mis)", 
         Approach = "Subtraction", Reliability = .)

rel_ccfe_d_regr_subt

Regression of Congruency (composite effect) in d [regression of regression]

rel_ccfe_d_regr_regr <- rel_reg(rel_ccf_d_cong_a_regr, rel_ccf_d_cong_m_regr, 
                                ccfe_d$cong_ali_regr, ccfe_d$cong_mis_regr) %>% 
  tibble(DV = "d", Condition = "(con_ali ~ inc_ali) ~ (con_mis ~ inc_mis)", 
         Approach = "Regression", Reliability = .)

rel_ccfe_d_regr_regr

Subtraction of RT [subtraction of subtraction]

rel_ccfe_rt_subt_subt <- rel_sub(rel_ccf_rt_cong_a_subt, rel_ccf_rt_cong_m_subt,
                                 ccfe_rt$cong_ali_subt, ccfe_rt$cong_mis_subt) %>% 
  tibble(DV = "rt", Condition = "(con_ali - inc_ali) - (con_mis - inc_mis)", 
         Approach = "Subtraction", Reliability = .)

rel_ccfe_rt_subt_subt

Regression of RT [regression of subtraction]

rel_ccfe_rt_subt_regr <- rel_reg(rel_ccf_rt_cong_a_subt, rel_ccf_rt_cong_m_subt,
                                 ccfe_rt$cong_ali_subt, ccfe_rt$cong_mis_subt) %>% 
  tibble(DV = "rt", Condition = "(con_ali - inc_ali) ~ (con_mis - inc_mis)", 
         Approach = "Regression", Reliability = .)

rel_ccfe_rt_subt_regr

Regression of RT [subtraction of regression]

rel_ccfe_rt_regr_subt <- rel_sub(rel_ccf_rt_cong_a_regr, rel_ccf_rt_cong_m_regr,
                                 ccfe_rt$cong_ali_regr, ccfe_rt$cong_mis_regr) %>% 
  tibble(DV = "rt", Condition = "(con_ali ~ inc_ali) - (con_mis ~ inc_mis)", 
         Approach = "Subtraction", Reliability = .)

rel_ccfe_rt_regr_subt

Regression of RT [regression of regression]

rel_ccfe_rt_regr_regr <- rel_reg(rel_ccf_rt_cong_a_regr, rel_ccf_rt_cong_m_regr,
                                 ccfe_rt$cong_ali_regr, ccfe_rt$cong_mis_regr) %>% 
  tibble(DV = "rt", Condition = "(con_ali ~ inc_ali) ~ (con_mis ~ inc_mis)", 
         Approach = "Regression", Reliability = .)

rel_ccfe_rt_regr_regr

3.2.3.7 Interim summary

rel_sum_ccf <- 
  bind_rows(rel_ccf_beha_ca, rel_ccf_beha_ia, rel_ccf_beha_cm, rel_ccf_beha_im, rel_ccf_beha_iso, # each condition
            rel_ccf_rt_ca, rel_ccf_rt_ia, rel_ccf_rt_cm, rel_ccf_rt_im, rel_ccf_rt_iso, # each condition
            rel_ccf_d_ca, rel_ccf_d_ia, rel_ccf_d_cm, rel_ccf_d_im, rel_ccf_d_iso, # each condition
            rel_fac_d_iso_subt, rel_fac_d_iso_regr, rel_fac_rt_iso_subt, rel_fac_rt_iso_regr, # facilitation (isolated)
            rel_fac_d_subt, rel_fac_d_regr, rel_fac_rt_subt, rel_fac_rt_regr, # facilitation 
            rel_int_d_iso_subt, rel_int_d_iso_regr, rel_int_rt_iso_subt, rel_int_rt_iso_regr, # interference (isolated)
            rel_int_d_subt, rel_int_d_regr, rel_int_rt_subt, rel_int_rt_regr, # interference 
            rel_ccf_d_cong_a_subt, rel_ccf_d_cong_m_subt, rel_ccf_rt_cong_a_subt, rel_ccf_rt_cong_m_subt, # congruency effects (subtraction)
            rel_ccf_d_cong_a_regr, rel_ccf_d_cong_m_regr, rel_ccf_rt_cong_a_regr, rel_ccf_rt_cong_m_regr, # congruency effects (subtraction)
            rel_ccfe_d_subt_subt, rel_ccfe_d_regr_subt, rel_ccfe_d_subt_regr, rel_ccfe_d_regr_regr, # composite effects
            rel_ccfe_rt_subt_subt, rel_ccfe_rt_regr_subt, rel_ccfe_rt_subt_regr, rel_ccfe_rt_regr_regr # composite effects
  ) %>% 
  mutate(Task = "CCF", .before=1)
rel_sum_ccf

3.3 Correlations among holistic processing effects

# custom function to print the two-sided correlation results via library(psych)
psych_cor2 <- function(.data, ...) {
  
  # .data: data frame only includes columns for correlation (no subject code column)
  # Usage: psych_cor2(hpe_acc_subt)
  
  # run correlations
  cor_tmp <- psych::corr.test(.data, ...)
  
  cor_tmp$ci %>% 
    rename(`lower95`=lower,
           `upper95`=upper) %>% 
    mutate(N = cor_tmp$n,
           p.adj = cor_tmp$p.adj,
           adjust = cor_tmp$adjust,
           alt = "two.sided") %>% 
    bind_cols(cor_tmp$ci.adj %>% mutate(r=1) %>% select(-r)) %>% 
    return()
}

3.3.1 Subtraction

3.3.1.1 Accuracy & d’

hpe_acc_subt <- pwe_acc %>% 
  inner_join(scfe_acc, by="Subject") %>% 
  inner_join(ccfe_d, by="Subject") %>% 
  transmute(pwe_subt, scfe_subt=-scfe_subt, ccfe_subt_subt) 

hpe_rt_subt <- pwe_rt %>% 
  inner_join(scfe_rt, by="Subject") %>% 
  inner_join(ccfe_rt, by="Subject") %>% 
  transmute(pwe_subt, scfe_subt=-scfe_subt, ccfe_subt_subt)

Correlations among holistic effects (accuracy or d’ with subtraction):

corPlot(hpe_acc_subt)

Correlations among holistic processing effects of accuracy (subtraction):

cor_hpe_acc_subt <- psych_cor2(hpe_acc_subt) %>% 
  rownames_to_column("pair") %>% 
  mutate(DV = "acc", 
         pair = c("pwe_subt ~ scfe_subt", "ccfe_subt ~ pwe_subt", "ccfe_subt ~ scfe_subt"),
         up_bound = c(
           upper_boundary(rel_pwe_acc_subt, rel_scfe_acc_subt),
           upper_boundary(rel_pwe_acc_subt, rel_ccfe_d_subt_subt),
           upper_boundary(rel_scfe_acc_subt, rel_ccfe_d_subt_subt)
         ), .before = 1) 
  
cor_hpe_acc_subt
plot_pwe_scfe_acc_subt <- ggplot(hpe_acc_subt, aes(x=pwe_subt, y=scfe_subt)) +
  geom_point(color=dark_colors[3], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Part-whole difference scores (accuracy)", 
       y = "Standard composite task difference scores (accuracy)") +
  coord_cartesian(ylim = c(-.35, .6)) +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")",
                          cor_hpe_acc_subt$r[1], cor_hpe_acc_subt$p[1]))
plot_pwe_scfe_acc_subt

plot_ccfe_pwe_acc_subt <- ggplot(hpe_acc_subt, aes(x=ccfe_subt_subt, y=pwe_subt)) +
  geom_point(color=dark_colors[1], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Complete composite task difference scores (d')", 
       y = "Part-whole difference scores (accuracy)") +
  # coord_cartesian(ylim = c(-.35, .6)) +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          cor_hpe_acc_subt$r[2], cor_hpe_acc_subt$p[2]))
plot_ccfe_pwe_acc_subt

plot_ccfe_scfe_acc_subt <- ggplot(hpe_acc_subt, aes(x=ccfe_subt_subt, y=scfe_subt)) +
  geom_point(color=dark_colors[2], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Complete composite task difference scores (d')", 
       y = "Standard composite task difference scores (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.3f\")", 
                          cor_hpe_acc_subt$r[3], cor_hpe_acc_subt$p[3]))
plot_ccfe_scfe_acc_subt

3.3.1.2 RT

Correlations among holistic effects (RT):

corPlot(hpe_rt_subt)

Correlations among holistic processing effects of RT (subtraction)

cor_hpe_rt_subt <- psych_cor2(hpe_rt_subt) %>% 
  rownames_to_column("pair") %>% 
  mutate(DV = "rt", 
         pair = c("pwe_subt ~ scfe_subt", "ccfe_subt ~ pwe_subt", "ccfe_subt ~ scfe_subt"),
         up_bound = c(
           upper_boundary(rel_pwe_rt_subt, rel_scfe_rt_subt),
           upper_boundary(rel_pwe_rt_subt, rel_ccfe_rt_subt_subt),
           upper_boundary(rel_scfe_rt_subt, rel_ccfe_rt_subt_subt)
         ), .before = 1) 
cor_hpe_rt_subt
plot_pwe_scfe_rt_subt <- ggplot(hpe_rt_subt, aes(x=pwe_subt, y=scfe_subt)) +
  geom_point(color=dark_colors[3], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Part-whole difference scores (RT)", 
       y = "Standard composite task difference scores (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          cor_hpe_rt_subt$r[1], cor_hpe_rt_subt$p[1]))

plot_pwe_scfe_rt_subt

plot_ccfe_pwe_rt_subt <- ggplot(hpe_rt_subt, aes(x=ccfe_subt_subt, y=pwe_subt)) +
  geom_point(color=dark_colors[1], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Complete composite task difference scores (RT)", 
       y = "Part-whole difference scores (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          cor_hpe_rt_subt$r[2], cor_hpe_rt_subt$p[2]))
plot_ccfe_pwe_rt_subt

plot_ccfe_scfe_rt_subt <- ggplot(hpe_rt_subt, aes(x=ccfe_subt_subt, y=scfe_subt)) +
  geom_point(color=dark_colors[2], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Complete composite task difference scores (RT)", 
       y = "Standard composite task difference scores (RT)") +
  coord_cartesian(ylim = c(-600, 1000)) +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          cor_hpe_rt_subt$r[3], cor_hpe_rt_subt$p[3]))
plot_ccfe_scfe_rt_subt

3.3.1.3 Plot

ggarrange(plot_pwe_scfe_acc_subt, plot_ccfe_pwe_acc_subt, plot_ccfe_scfe_acc_subt,
          plot_pwe_scfe_rt_subt, plot_ccfe_pwe_rt_subt, plot_ccfe_scfe_rt_subt,
          nrow = 2, ncol = 3)

3.3.2 Regression

3.3.2.1 Accuracy & d’

hpe_acc_regr <- pwe_acc %>% 
  inner_join(scfe_acc, by="Subject") %>% 
  inner_join(ccfe_d, by="Subject") %>% 
  transmute(pwe_regr, scfe_regr=-scfe_regr, ccfe_subt_regr)

hpe_rt_regr <- pwe_rt %>% 
  inner_join(scfe_rt, by="Subject") %>% 
  inner_join(ccfe_rt, by="Subject") %>% 
  transmute(pwe_regr, scfe_regr=-scfe_regr, ccfe_subt_regr)

Correlations among holistic effects (accuracy or d’ with regression):

corPlot(hpe_acc_regr)

Correlations among holistic processing effects of accuracy (regression)

cor_hpe_acc_regr <- psych_cor2(hpe_acc_regr) %>% 
  rownames_to_column("pair") %>% 
  mutate(DV = "acc", 
         pair = c("pwe_regr ~ scfe_regr", "ccfe_regr ~ pwe_regr", "ccfe_regr ~ scfe_regr"),
         up_bound = c(
           upper_boundary(rel_pwe_acc_regr, rel_scfe_acc_regr),
           upper_boundary(rel_pwe_acc_regr, rel_ccfe_d_subt_regr),
           upper_boundary(rel_scfe_acc_regr, rel_ccfe_d_subt_regr)
         ), .before = 1) 
cor_hpe_acc_regr
plot_pwe_scfe_acc_regr <- ggplot(hpe_acc_regr, aes(x=pwe_regr, y=scfe_regr)) +
  geom_point(color=dark_colors[3], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Part-whole residuals (accuracy)", 
       y = "Standard composite task residuals (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          cor_hpe_acc_regr$r[1], cor_hpe_acc_regr$p[1]))
plot_pwe_scfe_acc_regr

plot_ccfe_pwe_acc_regr <- ggplot(hpe_acc_regr, aes(x=ccfe_subt_regr, y=pwe_regr)) +
  geom_point(color=dark_colors[1], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Complete composite task residuals (d')", 
       y = "Part-whole residuals (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.3f\")", 
                          cor_hpe_acc_regr$r[2], cor_hpe_acc_regr$p[2]))
plot_ccfe_pwe_acc_regr

plot_ccfe_scfe_acc_regr <- ggplot(hpe_acc_regr, aes(x=ccfe_subt_regr, y=scfe_regr)) +
  geom_point(color=dark_colors[2], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Complete composite task residuals (d')", 
       y = "Standard composite task residuals (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.4f\")", 
                          cor_hpe_acc_regr$r[3], cor_hpe_acc_regr$p[3]))
plot_ccfe_scfe_acc_regr

3.3.2.2 RT

Correlations among holistic effects (RT with regression):

corPlot(hpe_rt_regr)

Correlations among holistic processing effects of RT (regression)

cor_hpe_rt_regr <- psych_cor2(hpe_rt_regr) %>% 
  rownames_to_column("pair") %>% 
  mutate(DV = "rt", 
         pair = c("pwe_regr ~ scfe_regr", "ccfe_regr ~ pwe_regr", "ccfe_regr ~ scfe_regr"),
         up_bound = c(
           upper_boundary(rel_pwe_rt_regr, rel_scfe_rt_regr),
           upper_boundary(rel_pwe_rt_regr, rel_ccfe_rt_subt_regr),
           upper_boundary(rel_scfe_rt_regr, rel_ccfe_rt_subt_regr)
         ), .before = 1) 
cor_hpe_rt_regr
plot_pwe_scfe_rt_regr <- ggplot(hpe_rt_regr, aes(x=pwe_regr, y=scfe_regr)) +
  geom_point(color=dark_colors[3], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Part-whole residuals (RT)", 
       y = "Standard composite task residuals (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          cor_hpe_rt_regr$r[1], cor_hpe_rt_regr$p[1]))
plot_pwe_scfe_rt_regr

plot_ccfe_pwe_rt_regr <- ggplot(hpe_rt_regr, aes(x=ccfe_subt_regr, y=pwe_regr)) +
  geom_point(color=dark_colors[1], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Complete composite task residuals (RT)", 
       y = "Part-whole residuals (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          cor_hpe_rt_regr$r[2], cor_hpe_rt_regr$p[2]))
plot_ccfe_pwe_rt_regr

plot_ccfe_scfe_rt_regr <- ggplot(hpe_rt_regr, aes(x=ccfe_subt_regr, y=scfe_regr)) +
  geom_point(color=dark_colors[2], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Complete composite task residuals (RT)", 
       y = "Standard composite task residuals (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          cor_hpe_rt_regr$r[3], cor_hpe_rt_regr$p[3]))
plot_ccfe_scfe_rt_regr

3.3.2.3 Plot

ggarrange(plot_pwe_scfe_acc_regr, plot_ccfe_pwe_acc_regr, plot_ccfe_scfe_acc_regr,
          plot_pwe_scfe_rt_regr, plot_ccfe_pwe_rt_regr, plot_ccfe_scfe_rt_regr,
          nrow = 2, ncol = 3)

3.3.3 Interim summary

cor_hpe <- bind_rows(cor_hpe_acc_subt, cor_hpe_rt_subt, 
                     cor_hpe_acc_regr, cor_hpe_rt_regr)
cor_hpe

3.4 Correlations for facilitation and interference

# custom function to print the two-sided correlation results via library(psych)
psych_cor1 <- function(x, y, ...) {
  
  # .data: data frame only includes columns for correlation (no subject code column)
  # Usage: psych_cor1(hpe_acc_subt)
  
  # run correlations
  cor_tmp <- psych::corr.test(x, y, alpha=.1, ...) # one-sided test
  
  rownames(cor_tmp) <- NULL
  
  cor_tmp$ci %>% 
    rename(`lower90`=lower,
           `upper90`=upper) %>% 
    mutate(N = cor_tmp$n,
           p = p/2, # one-sided p-value
           alt = "greater") %>% 
    # bind_cols(cor_tmp$ci.adj %>% mutate(r=1) %>% select(-r)) %>% 
    return()
}
# calculate the adjusted CI
psych_cor1_adj_CI <- function(.data, alpha=.10){
  # codes are modified from psych::corr.test()
  z <- fisherz(.data$r)
  ord <- order(z, decreasing = FALSE)
  
  dif.corrected <- qnorm(1-alpha/(2*(order(ord))))
  se <- 1/sqrt(.data$N - 3)
  
  # dif <- qnorm(1-alpha/2)
  # lower <- fisherz2r(z - dif * se)
  # upper <- fisherz2r(z + dif * se)
  .data %>% 
    mutate(lower.adj=fisherz2r(z - dif.corrected * se),
           upper.adj=fisherz2r(z + dif.corrected * se)) %>% 
    return()
}
hpe_acc <- pwe_acc %>% 
  inner_join(scfe_acc, by="Subject") %>% 
  inner_join(ccfe_d, by="Subject")

hpe_rt <- pwe_rt %>% 
  left_join(scfe_rt, by="Subject") %>% 
  left_join(ccfe_rt, by="Subject")

3.4.1 PW and Facilitation in CCF

3.4.1.1 Subtraction

accuracy

r_pwe_fac_acc_subt <- psych_cor1(pwe_acc$pwe_subt, ccfe_d$fac_subt) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "pwe_subt ~ fac_subt",
         DV = "acc",
         up_bound = upper_boundary(rel_pwe_acc_subt, rel_fac_d_subt),
         .before = 1)
r_pwe_fac_acc_subt
ggplot(hpe_acc, aes(x=fac_subt, y=pwe_subt)) +
  geom_point(color=dark_colors[1], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation difference scores (d')",
       y = "Part-whole difference scores (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f \")",
                          r_pwe_fac_acc_subt$r, r_pwe_fac_acc_subt$p))

accuracy isolated

r_pwe_fac_iso_acc_subt <- psych_cor1(pwe_acc$pwe_subt, ccfe_d$fac_subt_iso) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "pwe_subt ~ fac_subt_iso",
         DV = "acc",
         up_bound = upper_boundary(rel_pwe_acc_subt, rel_fac_d_iso_subt),
         .before = 1)
r_pwe_fac_iso_acc_subt
ggplot(hpe_acc, aes(x=fac_subt_iso, y=pwe_subt)) +
  geom_point(color=dark_colors[1], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation difference scores (d')",
       y = "Part-whole difference scores (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          r_pwe_fac_iso_acc_subt$r, 
                          r_pwe_fac_iso_acc_subt$p))

RT

r_pwe_fac_rt_subt <- psych_cor1(pwe_rt$pwe_subt, ccfe_rt$fac_subt) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "pwe_subt ~ fac_subt",
         DV = "rt", 
         up_bound = upper_boundary(rel_pwe_rt_subt, rel_fac_rt_subt),
         .before = 1)
r_pwe_fac_rt_subt
ggplot(hpe_rt, aes(x=fac_subt, y=pwe_subt)) +
  geom_point(color=dark_colors[1], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation difference scores (RT)",
       y = "Part-whole difference scores (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          r_pwe_fac_rt_subt$r, 
                          r_pwe_fac_rt_subt$p))

RT isolated

r_pwe_fac_iso_rt_subt <- psych_cor1(pwe_rt$pwe_subt, ccfe_rt$fac_subt_iso) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "pwe_subt ~ fac_subt_iso",
         DV = "rt", 
         up_bound = upper_boundary(rel_pwe_rt_subt, rel_fac_rt_iso_subt),
         .before = 1)
r_pwe_fac_iso_rt_subt
ggplot(hpe_rt, aes(x=fac_subt_iso, y=pwe_subt)) +
  geom_point(color=dark_colors[1], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation difference scores (RT)",
       y = "Part-whole difference scores (RT)") +
  coord_cartesian(ylim = c(-750, 1200)) +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          r_pwe_fac_iso_rt_subt$r, 
                          r_pwe_fac_iso_rt_subt$p))

3.4.1.2 Regressions

Accuracy

r_pwe_fac_acc_regr <- psych_cor1(pwe_acc$pwe_regr, ccfe_d$fac_regr) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "pwe_regr ~ fac_regr",
         DV = "acc", 
         up_bound = upper_boundary(rel_pwe_acc_regr, rel_fac_d_regr),
         .before = 1)
r_pwe_fac_acc_regr
ggplot(hpe_acc, aes(x=fac_regr, y=pwe_regr)) +
  geom_point(color=dark_colors[1], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation residuals (d')",
       y = "Part-whole residuals (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          r_pwe_fac_acc_regr$r, 
                          r_pwe_fac_acc_regr$p))

Accuracy isolated

r_pwe_fac_iso_acc_regr <- psych_cor1(pwe_acc$pwe_regr, ccfe_d$fac_regr_iso) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "pwe_regr ~ fac_regr_iso",
         DV = "acc", 
         up_bound = upper_boundary(rel_pwe_acc_regr, rel_fac_d_iso_regr),
         .before = 1)
r_pwe_fac_iso_acc_regr
ggplot(hpe_acc, aes(x=fac_regr_iso, y=pwe_regr)) +
  geom_point(color=dark_colors[1], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation residuals (d')",
       y = "Part-whole residuals (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          r_pwe_fac_iso_acc_regr$r, 
                          r_pwe_fac_iso_acc_regr$p))

RT

r_pwe_fac_rt_regr <- psych_cor1(pwe_rt$pwe_regr, ccfe_rt$fac_regr) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "pwe_regr ~ fac_regr",
         DV = "rt", 
         up_bound = upper_boundary(rel_pwe_rt_regr, rel_fac_rt_regr),
         .before = 1)
r_pwe_fac_rt_regr
ggplot(hpe_rt, aes(x=fac_regr, y=pwe_regr)) +
  geom_point(color=dark_colors[1], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation residuals (RT)",
       y = "Part-whole residuals (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          r_pwe_fac_rt_regr$r, 
                          r_pwe_fac_rt_regr$p))

RT isolated

r_pwe_fac_iso_rt_regr <- psych_cor1(pwe_rt$pwe_regr, ccfe_rt$fac_regr_iso) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "pwe_regr ~ fac_regr_iso",
         DV = "rt", 
         up_bound = upper_boundary(rel_pwe_rt_regr, rel_fac_rt_iso_regr),
         .before = 1)
r_pwe_fac_iso_rt_regr
ggplot(hpe_rt, aes(x=fac_regr_iso, y=pwe_regr)) +
  geom_point(color=dark_colors[1], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation residuals (RT)",
       y = "Part-whole residuals (RT)") +
  coord_cartesian(ylim = c(-650, 1200)) +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          r_pwe_fac_iso_rt_regr$r, 
                          r_pwe_fac_iso_rt_regr$p))

3.4.1.3 Interim summary

r_sum_pwe_fac <- bind_rows(r_pwe_fac_acc_subt, r_pwe_fac_iso_acc_subt,
                           r_pwe_fac_rt_subt, r_pwe_fac_iso_rt_subt,
                           r_pwe_fac_acc_regr, r_pwe_fac_iso_acc_regr,
                           r_pwe_fac_rt_regr, r_pwe_fac_iso_rt_regr)

Corrections for the baseline of isolated:

r_pwe_fac_iso <- r_sum_pwe_fac %>% 
  filter(endsWith(pair, "_iso")) %>% 
  mutate(p.adj = p.adjust(p, method = "holm"),
         adjust="holm") %>%  # apply Holm-Bonferroni
  psych_cor1_adj_CI() %>% 
  mutate(effect = "Facilitation")
r_pwe_fac_iso
plot_pwe_fac_acc_subt_iso <- ggplot(hpe_acc, aes(x=fac_subt_iso, y=pwe_subt)) +
  geom_point(color=dark_colors[1], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation difference scores (d')",
       y = "Part-whole difference scores (accuracy)") +
  coord_cartesian(ylim = c(-.25, .36)) +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.2f, \"
                          ,italic(p)[adj], \" = %0.2f\")", 
                          r_pwe_fac_iso_acc_subt$r, 
                          r_pwe_fac_iso_acc_subt$p,
                          r_pwe_fac_iso$p.adj[1]))
plot_pwe_fac_acc_subt_iso

plot_pwe_fac_rt_subt_iso <- ggplot(hpe_rt, aes(x=fac_subt_iso, y=pwe_subt)) +
  geom_point(color=dark_colors[1], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation difference scores (RT)",
       y = "Part-whole difference scores (RT)") +
  coord_cartesian(ylim = c(-770, 1200)) +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.2f, \"
                          ,italic(p)[adj], \" = %0.2f\")", 
                          r_pwe_fac_iso_rt_subt$r, 
                          r_pwe_fac_iso_rt_subt$p,
                          r_pwe_fac_iso$p.adj[2]))
plot_pwe_fac_rt_subt_iso

plot_pwe_fac_acc_regr_iso <-ggplot(hpe_acc, aes(x=fac_regr_iso, y=pwe_regr)) +
  geom_point(color=dark_colors[1], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation residuals (d')",
       y = "Part-whole residuals (accuracy)") +
  coord_cartesian(ylim = c(-.35, .19)) +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.2f, \"
                          ,italic(p)[adj], \" = %0.2f\")", 
                          r_pwe_fac_iso_acc_regr$r, 
                          r_pwe_fac_iso_acc_regr$p,
                          r_pwe_fac_iso$p.adj[3]))
plot_pwe_fac_acc_regr_iso

plot_pwe_fac_rt_regr_iso <- ggplot(hpe_rt, aes(x=fac_regr_iso, y=pwe_regr)) +
  geom_point(color=dark_colors[1], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation residuals (RT)",
       y = "Part-whole residuals (RT)") +
  coord_cartesian(ylim = c(-650, 1200)) +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.2f, \"
                          ,italic(p)[adj], \" = %0.2f\")", 
                          r_pwe_fac_iso_rt_regr$r, 
                          r_pwe_fac_iso_rt_regr$p,
                          r_pwe_fac_iso$p.adj[4]))
plot_pwe_fac_rt_regr_iso

ggarrange(plot_pwe_fac_acc_subt_iso, plot_pwe_fac_rt_subt_iso,
          plot_pwe_fac_acc_regr_iso, plot_pwe_fac_rt_regr_iso,
          nrow=2, ncol=2)

Corrections for the baseline of misaligned:

r_pwe_fac <- r_sum_pwe_fac %>% 
  filter(!endsWith(pair, "_iso")) %>% 
  mutate(p.adj = p.adjust(p, method = "holm"),
         adjust="holm") %>%  # apply Holm-Bonferroni
  psych_cor1_adj_CI()
r_pwe_fac
plot_pwe_fac_acc_subt <- ggplot(hpe_acc, aes(x=fac_subt, y=pwe_subt)) +
  geom_point(color=dark_colors[1], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation difference scores (d')",
       y = "Part-whole difference scores (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.2f, \"
                          ,italic(p)[adj], \" = %0.2f\")", 
                          r_pwe_fac_acc_subt$r, 
                          r_pwe_fac_acc_subt$p,
                          r_pwe_fac$p.adj[1]))
plot_pwe_fac_acc_subt

plot_pwe_fac_rt_subt <- ggplot(hpe_rt, aes(x=fac_subt, y=pwe_subt)) +
  geom_point(color=dark_colors[1], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation difference scores (RT)",
       y = "Part-whole difference scores (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.2f, \"
                          ,italic(p)[adj], \" = %0.2f\")",
                          r_pwe_fac_rt_subt$r, 
                          r_pwe_fac_rt_subt$p,
                          r_pwe_fac$p.adj[2]))
plot_pwe_fac_rt_subt

plot_pwe_fac_acc_regr <- ggplot(hpe_acc, aes(x=fac_regr, y=pwe_regr)) +
  geom_point(color=dark_colors[1], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation residuals (RT)",
       y = "Part-whole residuals (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.2f, \"
                          ,italic(p)[adj], \" = %0.2f\")",
                          r_pwe_fac_acc_regr$r, 
                          r_pwe_fac_acc_regr$p,
                          r_pwe_fac$p.adj[3]))
plot_pwe_fac_acc_regr

plot_pwe_fac_rt_regr <- ggplot(hpe_rt, aes(x=fac_regr, y=pwe_regr)) +
  geom_point(color=dark_colors[1], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Facilitation residuals (RT)",
       y = "Part-whole residuals (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.2f, \"
                          ,italic(p)[adj], \" = %0.2f\")", 
                          r_pwe_fac_rt_regr$r, 
                          r_pwe_fac_rt_regr$p,
                          r_pwe_fac$p.adj[4]))
plot_pwe_fac_rt_regr

ggarrange(plot_pwe_fac_acc_subt, plot_pwe_fac_rt_subt,
          plot_pwe_fac_acc_regr, plot_pwe_fac_rt_regr,
          nrow=2, ncol=2)

3.4.2 SCF and Interference in CCF

3.4.2.1 Subtraction

Accuracy

r_scfe_int_acc_subt <- psych_cor1(scfe_acc$scfe_subt, ccfe_d$int_subt) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "scfe_subt ~ int_subt",
         DV = "acc",
         up_bound = upper_boundary(rel_scfe_acc_subt, rel_int_d_subt),
         .before = 1)
r_scfe_int_acc_subt
ggplot(hpe_acc, aes(x=int_subt, y=scfe_subt)) +
  geom_point(color=dark_colors[2], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference difference scores (d')",
       y = "Standard composite task difference scores (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.4f\")", 
                          r_scfe_int_acc_subt$r, 
                          r_scfe_int_acc_subt$p))

Accuracy isolated

r_scfe_int_iso_acc_subt <- psych_cor1(scfe_acc$scfe_subt, ccfe_d$int_subt_iso) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "scfe_subt ~ int_subt_iso",
         DV = "acc",
         up_bound = upper_boundary(rel_scfe_acc_subt, rel_int_d_iso_subt),
         .before = 1)
r_scfe_int_iso_acc_subt
ggplot(hpe_acc, aes(x=int_subt_iso, y=scfe_subt)) +
  geom_point(color=dark_colors[2], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference difference scores (d')",
       y = "Standard composite task difference scores (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.5f\")", 
                          r_scfe_int_iso_acc_subt$r, 
                          r_scfe_int_iso_acc_subt$p))

RT

r_scfe_int_rt_subt <- psych_cor1(scfe_rt$scfe_subt, ccfe_rt$int_subt) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "scfe_subt ~ int_subt",
         DV = "rt",
         up_bound = upper_boundary(rel_scfe_rt_subt, rel_int_rt_subt),
         .before = 1)
r_scfe_int_rt_subt
ggplot(hpe_rt, aes(x=int_subt, y=scfe_subt)) +
  geom_point(color=dark_colors[2], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference difference scores (d')",
       y = "Standard composite task difference scores (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          r_scfe_int_rt_subt$r, 
                          r_scfe_int_rt_subt$p))

RT isolated

r_scfe_int_iso_rt_subt <- psych_cor1(scfe_rt$scfe_subt, ccfe_rt$int_subt_iso) %>%
  rownames_to_column("pair") %>% 
  mutate(pair = "scfe_subt ~ int_subt_iso",
         DV = "rt",
         up_bound = upper_boundary(rel_scfe_rt_subt, rel_int_rt_iso_subt),
         .before = 1)
r_scfe_int_iso_rt_subt
ggplot(hpe_rt, aes(x=int_subt_iso, y=scfe_subt)) +
  geom_point(color=dark_colors[2], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference difference scores (RT)",
       y = "Standard composite task difference scores (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.3f\")",  
                          r_scfe_int_iso_rt_subt$r, 
                          r_scfe_int_iso_rt_subt$p))

3.4.2.2 Regressions

Accuracy

r_scfe_int_acc_regr <- psych_cor1(scfe_acc$scfe_regr, ccfe_d$int_regr) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "scfe_regr ~ int_regr",
         DV = "acc",
         up_bound = upper_boundary(rel_scfe_acc_regr, rel_int_d_regr),
         .before = 1)
r_scfe_int_acc_regr
ggplot(hpe_acc, aes(x=int_regr, y=scfe_regr)) +
  geom_point(color=dark_colors[2], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference residuals (d')",
       y = "Standard composite task residuals (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.5f\")", 
                          r_scfe_int_acc_regr$r, 
                          r_scfe_int_acc_regr$p))

Accuracy isolated

r_scfe_int_iso_acc_regr <- psych_cor1(scfe_acc$scfe_regr, ccfe_d$int_regr_iso) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "scfe_regr ~ int_regr_iso",
         DV = "acc",
         up_bound = upper_boundary(rel_scfe_acc_regr, rel_int_d_iso_regr),
         .before = 1)
r_scfe_int_iso_acc_regr
ggplot(hpe_acc, aes(x=int_regr_iso, y=scfe_regr)) +
  geom_point(color=dark_colors[2], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference residuals (d')",
       y = "Standard composite task residuals (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.5f\")", 
                          r_scfe_int_iso_acc_regr$r, 
                          r_scfe_int_iso_acc_regr$p))

RT

r_scfe_int_rt_regr <- psych_cor1(scfe_rt$scfe_regr, ccfe_rt$int_regr) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "scfe_regr ~ int_regr",
         DV = "rt",
         up_bound = upper_boundary(rel_scfe_rt_regr, rel_int_rt_regr),
         .before = 1)
r_scfe_int_rt_regr
ggplot(hpe_rt, aes(x=int_regr, y=scfe_regr)) +
  geom_point(color=dark_colors[2], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference residuals (RT)",
       y = "Standard composite task residuals (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          r_scfe_int_rt_regr$r, 
                          r_scfe_int_rt_regr$p))

RT isolated

r_scfe_int_iso_rt_regr <- psych_cor1(scfe_rt$scfe_regr, ccfe_rt$int_regr_iso) %>% 
  rownames_to_column("pair") %>% 
  mutate(pair = "scfe_regr ~ int_regr_iso",
         DV = "rt",
         up_bound = upper_boundary(rel_scfe_rt_regr, rel_int_rt_iso_regr),
         .before = 1)
r_scfe_int_iso_rt_regr
ggplot(hpe_rt, aes(x=int_regr, y=scfe_regr)) +
  geom_point(color=dark_colors[2], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference residuals (RT)",
       y = "Standard composite task residuals (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \",italic(p), \" = %0.2f\")", 
                          r_scfe_int_iso_rt_regr$r, 
                          r_scfe_int_iso_rt_regr$p))

3.4.2.3 Interim summary

r_sum_scfe_int <- bind_rows(r_scfe_int_acc_subt, r_scfe_int_iso_acc_subt,
                            r_scfe_int_rt_subt, r_scfe_int_iso_rt_subt,
                            r_scfe_int_acc_regr, r_scfe_int_iso_acc_regr,
                            r_scfe_int_rt_regr, r_scfe_int_iso_rt_regr) 

Corrections for the baseline of isolated:

r_scfe_int_iso <- r_sum_scfe_int %>% 
  filter(endsWith(pair, "_iso")) %>% 
  mutate(p.adj = p.adjust(p, method = "holm"),
         adjust="holm") %>%  # apply Holm-Bonferroni
  psych_cor1_adj_CI() %>% 
  mutate(effect = "Interference")
r_scfe_int_iso
plot_scfe_int_acc_subt_iso <- ggplot(hpe_acc, aes(x=int_subt_iso, y=scfe_subt)) +
  geom_point(color=dark_colors[2], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference difference scores (d')",
       y = "Standard composite task difference scores (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.5f, \"
                          ,italic(p)[adj], \" = %0.5f\")", 
                          r_scfe_int_iso_acc_subt$r, 
                          r_scfe_int_iso_acc_subt$p,
                          r_scfe_int_iso$p.adj[1]))
plot_scfe_int_acc_subt_iso

plot_scfe_int_rt_subt_iso <- ggplot(hpe_rt, aes(x=int_subt_iso, y=scfe_subt)) +
  geom_point(color=dark_colors[2], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference difference scores (RT)",
       y = "Standard composite task difference scores (RT)") +
  coord_cartesian(ylim = c(-1100, 600)) +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.3f, \"
                          ,italic(p)[adj], \" = %0.3f\")", 
                          r_scfe_int_iso_rt_subt$r, 
                          r_scfe_int_iso_rt_subt$p,
                          r_scfe_int_iso$p.adj[2]))
plot_scfe_int_rt_subt_iso

plot_scfe_int_acc_regr_iso <- ggplot(hpe_acc, aes(x=int_regr_iso, y=scfe_regr)) +
  geom_point(color=dark_colors[2], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference residuals (d')",
       y = "Standard composite task residuals (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.5f, \"
                          ,italic(p)[adj], \" = %0.5f\")", 
                          r_scfe_int_iso_acc_regr$r, 
                          r_scfe_int_iso_acc_regr$p,
                          r_scfe_int_iso$p.adj[3]))
plot_scfe_int_acc_regr_iso

plot_scfe_int_rt_regr_iso <- ggplot(hpe_rt, aes(x=int_regr, y=scfe_regr)) +
  geom_point(color=dark_colors[2], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference residuals (RT)",
       y = "Standard composite task residuals (RT)") +
  coord_cartesian(ylim = c(-730, 470)) +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.2f, \"
                          ,italic(p)[adj], \" = %0.2f\")", 
                          r_scfe_int_iso_rt_regr$r, 
                          r_scfe_int_iso_rt_regr$p,
                          r_scfe_int_iso$p.adj[4]))
plot_scfe_int_rt_regr_iso

ggarrange(plot_scfe_int_acc_subt_iso, plot_scfe_int_rt_subt_iso,
          plot_scfe_int_acc_regr_iso, plot_scfe_int_rt_regr_iso,
          nrow=2, ncol=2)

Corrections for the baseline of misaligned:

r_scfe_int <- r_sum_scfe_int %>% 
  filter(!endsWith(pair, "_iso")) %>% 
  mutate(p.adj = p.adjust(p, method = "holm"),
         adjust="holm") %>%  # apply Holm-Bonferroni
  psych_cor1_adj_CI()
r_scfe_int
plot_scfe_int_acc_subt <- ggplot(hpe_acc, aes(x=int_subt, y=scfe_subt)) +
  geom_point(color=dark_colors[2], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference difference scores (d')",
       y = "Standard composite task difference scores (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.4f, \"
                          ,italic(p)[adj], \" = %0.3f\")", 
                          r_scfe_int_acc_subt$r, 
                          r_scfe_int_acc_subt$p,
                          r_scfe_int$p.adj[1]))
plot_scfe_int_acc_subt

plot_scfe_int_rt_subt <- ggplot(hpe_rt, aes(x=int_subt, y=scfe_subt)) +
  geom_point(color=dark_colors[2], shape=subt_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference difference scores (d')",
       y = "Standard composite task difference scores (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.2f, \"
                          ,italic(p)[adj], \" = %0.2f\")", 
                          r_scfe_int_rt_subt$r, 
                          r_scfe_int_rt_subt$p,
                          r_scfe_int$p.adj[2]))
plot_scfe_int_rt_subt

plot_scfe_int_acc_regr <- ggplot(hpe_acc, aes(x=int_regr, y=scfe_regr)) +
  geom_point(color=dark_colors[2], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference residuals (d')",
       y = "Standard composite task residuals (accuracy)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.5f, \"
                          ,italic(p)[adj], \" = %0.5f\")", 
                          r_scfe_int_acc_regr$r, 
                          r_scfe_int_acc_regr$p,
                          r_scfe_int$p.adj[3]))
plot_scfe_int_acc_regr

plot_scfe_int_rt_regr <- ggplot(hpe_rt, aes(x=int_regr, y=scfe_regr)) +
  geom_point(color=dark_colors[2], shape=regr_shape) + 
  geom_smooth(color="black", method=lm, formula='y ~ x') +
  labs(x = "Interference residuals (RT)",
       y = "Standard composite task residuals (RT)") +
  annotate("text", x=Inf, y=-Inf, hjust = 1.1, vjust = -.5, parse = TRUE,
           label= sprintf("paste(italic(r), \" = %0.2f, \"
                          ,italic(p), \" = %0.2f, \"
                          ,italic(p)[adj], \" = %0.2f\")", 
                          r_scfe_int_rt_regr$r, 
                          r_scfe_int_rt_regr$p,
                          r_scfe_int$p.adj[4]))
plot_scfe_int_rt_regr

ggarrange(plot_scfe_int_acc_subt, plot_scfe_int_rt_subt,
          plot_scfe_int_acc_regr, plot_scfe_int_rt_regr,
          nrow=2, ncol=2)

3.5 All pairwise correlations (just for fun)

3.5.1 Behavioral choices

Correlations among effects (accuracy or d’ with subtraction):

pwe_acc %>% 
  inner_join(scfe_acc, by="Subject") %>% 
  inner_join(ccfe_d, by="Subject") %>% 
  select(-Subject, -contains("regr")) %>% 
  cor.plot()

Correlations among effects (accuracy or d’ with regression):

pwe_acc %>% 
  inner_join(scfe_acc, by="Subject") %>% 
  inner_join(ccfe_d, by="Subject") %>% 
  select(-Subject, -contains("subt"), ccfe_subt_regr) %>% 
  cor.plot()

3.5.2 Response times

Correlations among effects (RT with subtraction):

pwe_rt %>% 
  left_join(scfe_rt, by="Subject") %>% 
  left_join(ccfe_rt, by="Subject") %>% 
  select(-Subject, -contains("regr")) %>% 
  cor.plot()

Correlations among effects (RT with regression):

pwe_rt %>% 
  left_join(scfe_rt, by="Subject") %>% 
  left_join(ccfe_rt, by="Subject") %>% 
  select(-Subject, -contains("subt"), ccfe_subt_regr) %>% 
  cor.plot()

3.6 Publication

3.6.1 Table 1: reliability

rel_sum <- bind_rows(rel_sum_pw, rel_sum_scf, rel_sum_ccf) %>% 
  select(-niter)
# write_csv(rel_sum, "reliability_summary.csv")
rel_sum

3.6.2 Table 2: correlations

cor_sum <- cor_hpe %>% 
  rename(upper = upper95,
         lower = lower95) %>% 
  mutate(effect = "Holistic Processing", .before=1) %>%  
  mutate(CI = .95) %>% 
  bind_rows(bind_rows(r_pwe_fac_iso, r_scfe_int_iso) %>% 
              rename(upper = upper90,
                     lower = lower90) %>% 
              mutate(CI = .90))
# write_csv(cor_sum, "correlation_summary.csv")
cor_sum

3.6.3 Figure 3: correlations among holistic processing effects

plot_cor_hpe <- 
  ggarrange(plot_pwe_scfe_acc_subt, plot_ccfe_pwe_acc_subt, plot_ccfe_scfe_acc_subt,
            plot_pwe_scfe_rt_subt, plot_ccfe_pwe_rt_subt, plot_ccfe_scfe_rt_subt,
            plot_pwe_scfe_acc_regr, plot_ccfe_pwe_acc_regr, plot_ccfe_scfe_acc_regr,
            plot_pwe_scfe_rt_regr, plot_ccfe_pwe_rt_regr, plot_ccfe_scfe_rt_regr,
            nrow = 4, ncol = 3)

# plot_cor_hpe
plot_cor_hpe

3.6.4 Figure 4: correlations with facilitation and interference

plot_cor_fi <- ggarrange(plot_pwe_fac_acc_subt_iso, plot_pwe_fac_rt_subt_iso,
                         plot_pwe_fac_acc_regr_iso, plot_pwe_fac_rt_regr_iso,
                         plot_scfe_int_acc_subt_iso, plot_scfe_int_rt_subt_iso,
                         plot_scfe_int_acc_regr_iso, plot_scfe_int_rt_regr_iso,
                         nrow=4, ncol=2)
# plot_cor_fi
plot_cor_fi

4 Demographic information

# To run the code in this chunk, you need to put all demographic data files into data/.

# read demographic files
# df_demo <- list.files("data", pattern = "*.csv", full.names = TRUE) %>%
#   map_dfr(read_csv, col_types="cccccccccdidfcffffffff", .id="id") %>% 
#   select(ProlificID=`Participant id`, Age, Sex) %>% 
#   right_join(df_subjects, by="ProlificID") %>% 
#   select(Subject, Age, Sex)
# 
# saveRDS(df_demo, here("data", "df_demo.rds"))

Sex

df_demo <- readRDS(here("data", "df_demo.rds")) %>% 
  filter(!Subject %in% subjlist_ex) # exclude outlier subjects

df_demo %>% 
  mutate(isFemale = Sex=="Female",
         isMale = Sex=="Male") %>% 
  summarize(count_female = sum(isFemale),
            count_male = sum(isMale),
            count_all = n())

Age

df_demo %>% 
  summarize(Age_mean = mean(as.integer(Age)),
            Age_SD = sd(Age),
            Age_min = min(Age),
            Age_max = max(Age))

Session information

sessionInfo()
## R version 4.1.3 (2022-03-10)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur/Monterey 10.16
## 
## Matrix products: default
## BLAS:   /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.1/Resources/lib/libRlapack.dylib
## 
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] here_1.0.1       ggpubr_0.4.0     psych_2.2.5      emmeans_1.8.0    optimx_2022-4.30 lmerTest_3.1-3   lme4_1.1-30      Matrix_1.4-0     forcats_0.5.1    stringr_1.4.0    dplyr_1.0.9      purrr_0.3.4      readr_2.1.2      tidyr_1.2.0      tibble_3.1.8     ggplot2_3.3.6    tidyverse_1.3.2 
## 
## loaded via a namespace (and not attached):
##  [1] nlme_3.1-155        fs_1.5.2            lubridate_1.8.0     insight_0.18.2      httr_1.4.3          rprojroot_2.0.3     numDeriv_2016.8-1.1 tools_4.1.3         backports_1.4.1     bslib_0.4.0         utf8_1.2.2          R6_2.5.1            mgcv_1.8-39         DBI_1.1.3           colorspace_2.0-3    withr_2.5.0         gridExtra_2.3       tidyselect_1.1.2    mnormt_2.1.0        compiler_4.1.3      cli_3.3.0           rvest_1.0.2         xml2_1.3.3          labeling_0.4.2      bayestestR_0.12.1   sass_0.4.2          scales_1.2.0        mvtnorm_1.1-3       digest_0.6.29       minqa_1.2.4         rmarkdown_2.14      tinylabels_0.2.3    pkgconfig_2.0.3     htmltools_0.5.3     highr_0.9           dbplyr_2.2.1        fastmap_1.1.0       rlang_1.0.4         readxl_1.4.0        rstudioapi_0.13     farver_2.1.1        jquerylib_0.1.4     generics_0.1.3      jsonlite_1.8.0      car_3.1-0           googlesheets4_1.0.1 magrittr_2.0.3      parameters_0.18.2   Rcpp_1.0.9          munsell_0.5.0       fansi_1.0.3         abind_1.4-5         lifecycle_1.0.1     stringi_1.7.8       yaml_2.3.5          carData_3.0-5       MASS_7.3-55         grid_4.1.3          parallel_4.1.3      crayon_1.5.1        lattice_0.20-45     cowplot_1.1.1       haven_2.5.0         splines_4.1.3       hms_1.1.1           knitr_1.39          pillar_1.8.0        uuid_1.1-0          boot_1.3-28         estimability_1.4.1  ggsignif_0.6.3      effectsize_0.7.0.5  reprex_2.0.1        glue_1.6.2         
## [75] evaluate_0.16       modelr_0.1.8        vctrs_0.4.1         nloptr_2.0.3        tzdb_0.3.0          cellranger_1.1.0    papaja_0.1.1        gtable_0.3.0        datawizard_0.5.0    assertthat_0.2.1    cachem_1.0.6        xfun_0.32           xtable_1.8-4        broom_1.0.0         coda_0.19-4         rstatix_0.7.0       googledrive_2.0.0   gargle_1.2.0        ellipsis_0.3.2      xaringanExtra_0.7.0